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ABSTRACT 


A  theoretical  investigation  of  the  ducted  turbu¬ 
lent  mixing  of  coaxial  streams  is  presented.  The  effects  of 
equilibrium  chemical  reactions  in  the  mixing  layer  are  in¬ 
cluded  in  the  analysis.  In  Part  I,  the  theory  is  developed 
for  ducted  mixing  systems  in  which  the  inviscid  portions  of 
the  flow  are  assumed  to  be  one-dimensional.  The  integral 
method  is  used,  and  the  mixing  zone  velocity  profiles  are 
assumed  to  be  similar.  Turbulent  Prandtl  and  Lewis  numbers 
are  assumed  to  be  unity,  and  the  turbulent  shear  stress 
along  one  control  surface  in  the  mixing  layer  is  computed  by 
use  of  a  semi-empirical  model  for  the  eddy  viscosity.  The 
theory  is  formulated  to  include  the  downstream  regime  where 
the  mixing  layer  extends  across  the  entire  duct.  The  system 
of  integral  equations  is  transformed  into  a  form  suitable  for 
numerical  solution  by  standard  methods.  In  Part  II,  the  theory 
is  extended  to  include  the  influence  of  a  nonuniform  supersonic 
inviscid  core  flow.  A  technique  is  developed  for  coupling  the 
numerical  solution  of  the  inviscid  core  flow  (from  the  method 
of  characteristics)  with  the  solution  of  the  integral  equations 
for  the  mixing  layer.  The  theory  is  correlated  with  low  speed 
experiments,  with  supersonic  air-air  ejector  experiments  and 
with  experiments  on  ducted  rocket-air  mixing.  The  correlations 
show  that  the  theory  gives  satisfactory  results  for  the  per¬ 
formance  of  ducted  rocket-air  mixing  systems  which  are 
representative  of  air  augmented  rockets. 
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SECTION  1 
INTRODUCTION 

Ducted  turbulent  mixing  of  coaxial  streams  occurs 
in  many  devices  of  practical  interest.  Typical  examples  are 
the  jet  pump  (or  air— air  ejector)  and  composite  propulsion 
systems  such  as  the  air  augmented  rocket  (Refs.  1,  2).  For 
propulsive  applications  one  must  consider  the  possibility 
that  exothermic  chemical  reactions  will  occur  in  the  mixing 
layer . 

The  basic  objective  of  this  research  is  to  develop 
a  theory  to  describe  the  ducted  turbulent  mixing  of  coaxial 
streams,  Including  chemical  reactions.  Emphasis  has  been 
placed  on  relatively  long  mixing  systems  in  which  the  mixing 
layer  may  extend  over  most  or  all  of  the  duct  cross  section 
at  the  duct  exit  plane.  The  duct  pressure  distribution  will 
be  strongly  influenced  by  the  thick  mixing  layers,  and  will 
be  very  different  from  the  invlscid  pressure  distribution. 

In  other  words,  this  may  be  considered  a  "strong  viscous 
interaction"  problem. 

Consider  the  mixing  system  shown  schematically  in 
Fig.  1,  in  which  the  primary  and  secondary  fluids  are  speci¬ 
fied.  The  specific  objective  of  this  work  is  to  predict  the 
secondary  mass  flow  rate,  w^,  and  the  duct  wall  pressure 
distribution  if  the  following  parameters  are  prescribed: 

1 .  Geometry 

2.  Primary  (central)  stream  initial  conditions 

3.  Secondary  stream  stagnation  pressure,  p__ , 
and  stagnation  temperature,  T^^^ 

4.  Back  pressure,  p^ 

Three  distinct  flow  regimes  are  shown  in  the 
mixing  flow  field  of  Fig.  1.  In  the  first  regime,  turbulent 
mixing  occurs  between  the  inviscid  secondary  flow  and  the 


1 


N> 


Fig.  1  Schematic  of  Ducted  Mixing  System 
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core  of  invlscid  primary  flow.  In  the  second  regime,  the 
inviscid  core  has  been  dissipated,  but  a  region  of  invlscid 
secondary  flow  exists  near  the  duct  wall.  The  third  regime 
occurs  after  the  mixing  layer  has  spread  to  the  duct  wall, 
and  the  flow  is  entirely  turbulent.  The  mixing  Is  essen-> 
tially  free  turbulent  in  nature  in  the  first  and  second 
regimes,  although  not  Isobarlc. 

Several  modes  of  operation  are  possible  for  a 
ducted  mixing  system  when  the  primary  stream  Is  initially 
supersonic  and  the  secondary  stream  is  initially  subsonic. 
These  operating  modes  are  distinguished  by  the  factor  which 
limits  the  secondary  mass  flow  rate,  w^.  The  "upstream 
choking"  mode  occurs  when  the  back  pressure  is  relatively 
low,  and  the  secondary  stream  static  pressure  is  low 
relative  to  the  primary  flow  (Fig.  2a).  The  primary  stream 
expands,  causing  the  secondary  flow  to  choke  at  some  station 
near  the  inlet  of  the  mixing  duct,  and  the  secondary  flow 
is  supersonic  downstream  of  the  choking  station.  Fabri  and 
Paulon  (Ref.  3)  called  this  mode  the  "supersonic  regime”  in 
their  work  on  air-air  ejectors.  As  the  secondary  stream 
stagnation  pressure  is  increased,  the  choking  station  will 
move  upstream.  If  the  secondary  stagnation  pressure  is 
increased  sufficiently  the  choking  station  will  be  fixed  at 
the  initial  section,  with  the  primary  nozzle  either  correctly 
expanded  or  overexpanded.  Fabri  and  Paulon  called  this 
limiting  case  of  the  upstream  choking  mode  the  "saturated 
supersonic  regime."  The  secondary  stream  minimum  area  is 
a  function  only  of  the  geometry,  and  the  mass  flow,  w  ,  is 

SL 

directly  proportional  to  p^^  in  this  limiting  case. 

The  "downstream  choking"  mode  (Fig.  2b)  is  en¬ 
countered  when  the  back  pressure  is  low,  and  w  is  limited 

a 

by  choking  of  the  flow  at  or  near  the  duct  exit.  The 
secondary  flow  is  subsonic  throughout  the  duct.  This  mode 
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a.  Upstream  Choking  Mode 
Pig.  2  Operational  Modes  of  a  Ducted  Mixing  System 

t 


4 


AEDC-TR-68-270 


b.  Downstream  Choking  Mode 
Fig.  2  Continged 
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will  probably  not  occur  in  a  cylindrical  or  divergent  mixing 
duct  unless  the  primary  and  secondary  fluids  have  greatly 
different  densities,  or  unless  chemical  reactions  occur  in 
the  mixing  process.  The  downstream  choking  mode  is  impor¬ 
tant  for  certain  propulsive  applications. 

The  "back  pressure  dependent"  mode  of  operation 
(Fig.  2c)  occurs  when  the  back  pressure,  p^,  is  sufficiently 
high  to  unchoke  the  duct  flow.  The  secondary  flow  is  subsonic 
throughout  the  duct,  and  the  duct  exit  pressure  exactly 
matches  the  back  pressure.  Fabri  and  Paulon  called  this  mode 
the  "mixed  regime"  and  it  is  commonly  encountered  in  jet 
pump  applications. 

If  it  is  to  be  considered  adequate,  a  theory  for 
the  mixing  duct  flow  must  be  capable  of  predicting  the  per¬ 
formance  in  each  of  the  operational  modes  shown  in  Fig.  2. 

It  should  be  noted  that  certain  limiting  cases  of 
ducted  mixing  can  be  successfully  analyzed  by  application  of 
one— dimensional  theoretical  concepts  (Ref.  3,  for  example).- 
If  the  mixing  duct  is  relatively  long  and  of  constant  cross 
section,  and  the  system  is  operating  in  the  back  pressure 
dependent  mode,  then  the  performance  can  be  predicted  by 
assuming  that  the  flow  is  fully  mixed  and  uniform  at  the  duct 
exit.  If  the  initial  conditions  are  specified,  then  overall 
conservation  equations  for  momentum,  mass  flow  and  energy 
can  be  solved  for  the  one— dimensional  flow  conditions  at  the 
duct  exit.  This  one— dimensional  theory  gives  no  information 
about  the  duct  length  required  to  achieve  complete  mixing 
because  the  equations  include  no  characteristic  length 
scale.  Also,  the  one— dimensional  approach  is  not  applicable 
to  mixing  ducts  of  variable  cross  section  because  the  unknown 
axial  wall  pressure  forces  cannot  be  predicted.  The  detailed 
mixing  duct  flow  must  be  treated  in  order  to  compute  the  un¬ 
known  wall  pressure  integral  in  variable  area  configurations. 


7 


AEDC-TR-66-270 


Some  limited  success  has  been  achieved  in  analyzing 
constant  area  ejector  systems  operating  in  the  upstream 
choking  mode  by  application  of  "quasi-one-dimens ional"  con¬ 
cepts.  Fabri  and  Paulon  (Ref.  3)  and  Kochendorfer  and 
Rousso  (Ref.  4)  have  used  such  quasi— one-dimensional  approaches. 
Viscous  effects  are  neglected,  and  again,  no  characteristic 
length  appears  in  the  equations.  The  analysis  is  limited 
to  mixing  ducts  of  constant  cross  sectional  area.  The 
primary  stream  is  assumed  to  expand  isentropically  and  one— 
dimensionally,  while  the  secondary  stream  contracts.  Fabri 
and  Paulon  assumed  that  the  outer  stream  is  isentropic  and 
that  M  =  1  at  the  choking  section.  Solution  of  the  com— 
bined  momentum  equation  for  both  streams  makes  it  unnecessary 
to  assume  that  the  primary  and  secondary  pressures  are  equal 
at  the  choking  section.  Kochendorfer  and  Rousso  formulated 
their  analysis  somewhat  differently.  They  assumed  that 
=  1  and  that  the  primary  and  secondary  pressures  are 
equal  at  the  choking  section.  The  combined  stream  momentum 
equation  was  used  to  define  the  secondary  stream  total 
pressure  at  the  choking  section. 

The  key  assumption  in  these  quasi— one— dimensional 
approaches  is  that  the  area  of  the  primary  stream  at  the 
choking  section  is  the  same  as  the  one— dimensional  isentropic 
area  at  the  same  static  pressure.  This  assumption  is  approxi¬ 
mately  true  when  the  primary  nozzle  flow  is  nearly  uniform 
and  not  too  highly  under expanded,  but  can  be  seriously  in 
error  in  many  practical  cases.  The  neglect  of  viscous 
effects  is  less  serious  for  supersonic  ejectors  with  rela¬ 
tively  large  secondary  flows.  For  this  case,  Chow  and  Addy 
(Ref.  5)  have  shown  that  viscous  effects  are  small.  It 
should  be  noted,  however,  that  viscous  effects  are  important 
for  an  air— air  ejector  operating  with  low  or  zero  secondary 
flow  rates. 
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The  quasi— oae-diraensional  theories  are  simple  to 
apply  and  give  good  engineering  results  for  many  constant 
area  ejector  configurations.  These  theories,  however,  are 
too  limited  in  scope  for  application  to  many  ejector  problems. 

In  Part  I  of  this  report,  a  theoretical  model  is 
developed  for  the  ducted  mixing  system,  based  on  the  integral 
form  of  the  steady  boundary  layer  equations.  Both  the  in- 
viscid  primary  and  secondary  streams  are  assumed  to  be  one- 
dimensional.  The  mixing  zone  chemistry  is  assumed  to  be 
either  frozen  or  in  equilibrium.  The  turbulent  shear  along 
one  control  surface  in  the  mixing  layer  is  computed  by  means 
of  a  semi-empirical  function  for  the  eddy  viscosity.  This 
theory,  hereafter  called  the  "1— D  Core  Theory,”  agrees  well 
with  experimental  results  for  subsonic  mixing  flows,  but  is 
inadequate  when  the  primary  stream  is  supersonic  and  non- 
uniform. 

In  Part  II  of  this  report,  the  first  regime  theory 
is  extended  to  include  treatment  of  the  inviscid  core  flow 
with  the  irrotational  method  of  characteristics  (This  theory 
will  be  called  the  "2— D  Core  Theory.”).  Inclusion  of  the 
method  of  characteristics  greatly  complicates  the  solution; 
however,  the  modified  theory  can  be  applied  to  supersonic 
primary  flows  which  are  initially  nonuniform  (conical  pri¬ 
mary  nozzles,  for  example),  or  to  cases  where  the  primary 
nozzle  is  underexpanded. 

It  should  be  noted  that  the  2— D  Core  Theory  differs 
from  the  1— D  Core  Theory  only  in  the  first  regime  in  which 
there  is  an  inviscid  core  of  primary  flow  (Fig.  1).  The  two 
theories  are  identical  in  the  second  and  third  regimes, 
i.e.,  after  the  inviscid  core  has  been  dissipated. 

Most  of  the  comparisons  with  experimental  results 
have  been  made  with  the  2— D  Core  Theory.  The  1— D  Core  Theory 
has  been  correlated  only  with  some  low  speed  experimental 
results,  both  to  establish  the  empirical  incompressible  eddy 
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viscosity  constants  and  to  establish  the  validity  of  the 
theory  for  low  speed  flows. 

Many  assumptions  have  been  made  in  the  analytical 
model  developed  during  this  investigation.  Some  of  these 
assumptions  are  necessary  because  the  current  knowledge  of 
turbulent  flow  processes  is  very  meager,  even  for  flow  con¬ 
figurations  much  less  complex  than  the  one  considered  here. 
Other  assumptions  have  been  made  to  render  the  analytical 
model  tractable  from  a  mathematical  point  of  view.  Recog¬ 
nizing  the  many  simplifications  in  the  theoretical  model, 
it  is  essential  that  its  validity  be  established  by  detailed 
comparisons  with  experiment  in  the  flow  regime  of  interest. 
Extensive  experiments  have  been  made  by  the  author  and  his 
co-workers  at  AEDC  on  a  ducted  mixing  system  in  which  the 
primary  flow  was  the  exhaust  from  a  fuel-rich  rocket  (Refs. 

6,  7).  The  secondary  stream  was  room  temperature  air,  and 
various  mixing  duct  configurations  were  used  so  that  all  of 
the  operating  modes  shown  in  Fig.  2  were  encountered.  Cor¬ 
relations  of  the  2-D  Core  Theory  with  these  rocket-air  mixing 
experiments,  as  well  as  with  experimental  results  for  air- 
air  ejectors,  are  discussed  in  Sections  XVI  and  XVII. 
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SECTION  II 

SUMMARY  OF  ANALYTICAL  APPROACHES  FOR  TURBULENT  MIXING 

The  purpose  of  this  section  is  to  give  a  brief  sur¬ 
vey  of  the  various  methods  which  have  been  used  in  the  analysis 
of  free  turbulent  mixing,  and  to  discuss  the  methods  in  view 
of  the  current  state  of  experimental  knowledge.  A  complete 
review  of  the  extensive  mixing  literature  is  not  intended; 
the  objective  is  to  mention  the  various  analytical  methods 
so  that  the  approach  taken  in  this  investigation  can  be  viewed 
in  historical  perspective.  More  detailed  information  on  the 
earlier  work  in  this  field  can  be  found  in  the  excellent 
monograph  of  Abramovich  (Ref.  8)  and  in  Refs.  ^11. 

The  analysis  of  free  turbulent  mixing  is  neces¬ 
sarily  semi— empirical  because  no  fundamental  quantitative 
flow  model  for  turbulence  is  currently  available.  Numerous 
experimental  Investigations,  however,  have  established  cer¬ 
tain  approximations  which  provide,  reasonably  accurate  predic¬ 
tions  of  the  mean  flow  field,  at  least  for  relatively  simple 
flow  configurations.  All  of  the  analytical  approaches  have 
incorporated  the  simplification  that  the  complex  fluctuating 
turbulent  flow  field  can  be  replaced  with  an  equivalent  steady 
flow.  The  fluctuating  nature  of  the  flow  is  usually  treated 
by  use  of  "apparent”  turbulent  transport  properties  which 
are  analogous  to  molecular  transport  properties  in  laminar 
flow.  The  usual  boundary  layer  approximations  have  been 
used  in  most  analyses.  Most  of  the  existing  analyses  have 
been  made  for  isobarlc  flow,  l.e.,  Jets,  wakes  and  mixing 
layers  surrounded  by  uniform  inviscid  streams  of  large  ex¬ 
tent  . 

1.  Low  SpMd  Mixing  (Constant  Dens(ty) 

Models  for  turbulent  shear  —  In  1925  Prandtl 
formulated  his  famous  mixing  length  theory  for  turbulent 
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shear.  He  hypothesized  that  the  mean  value  of  the  fluc¬ 
tuating  velocity  component  in  the  turbulent  field  is  equal 
to  the  product  of  the  local  mean  velocity  gradient  times  a 
characteristic  mixing  length,  1.  The  mixing  length  is 
proportional  to  the  local  width  of  the  mixing  layer  for  free 
mixing  (.Ji  =  cb)  ,  and  is  assumed  to  be  constant  across  the 
layer.  A  small  element  of  fluid  is  assumed  to  conserve  its 
axial  velocity  while  moving  transversely  over  the  mixing 
length.  The  turbulent  shear  stress  (Reynolds  stress)  in 
the  axial  direction  is 

z-=  -  ririP  -  p J 

where  c  is  an  empirical  constant,  u  is  the  mean  axial 
velocity,  p  is  the  density,  b  is  the  width  of  the  mixing 
layer  and  u'  and  v'  are  the  local  components  of  the  fluc¬ 
tuating  velocity. 

Taylor’s  vorticity  transport  theory  is  based  on 
a  similar  mixing  length  concept,  with  the  vorticity  of  a 
fluid  element  being  conserved  over  the  mixing  length.  Both 
Prandtl's  and  Taylor's  theories  assume  that  the  eddy  length 
scale,  i,  is  much  less  than  the  local  width  of  the  mixing 
layer.  Taylor's  vorticity  transport  theory  results  in  the 
same  expression  for  turbulent  shear  as  Prandtl's  mixing 
length  theory  for  two-dimensional  mixing.  Because  of  its 
complexity  in  all  but  the  simplest  flows,  Taylor's  theory 
has  been  used  much  Less  frequently  than  the  Pcandtl  mixing 
length  theory. 

The  Boussinesq  concept  of  an  "apparent"  or  "eddy" 
viscosity  has  proven  to  be  very  useful  in  the  analysis  of 
free  mixing  problems.  The  axial  shear  stress  is 
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where  €  is  the  eddy  kinematic  viscosity.  The  eddy  viscosity 
is  not  a  fluid  property,  but  is  related  to  the  dynamics  of 
the  particular  problem. 

Prandtl's  mixing  length  theory  gives  the  fol¬ 
lowing  relation  for  the  eddy  viscosity: 


Thus,  Prandtl's  mixing  length  theory  states  that  the  eddy 
viscosity  varies  with  |du/Sy|  across  the  mixing  layer. 

In  1942  Prandtl  proposed  another  model  for  the 
turbulent  eddy  viscosity,  based  on  the  assumption  that  v'Jl 
(equivalent  to  e)  is  related  to  the  local  mean  velocity 
gradient  across  the  mixing  layer.  This  model,  which  implies 
that  the  turbulent  eddy  scale  is  of  the  same  order  as  the 
mixing  zone  width,  yields  the  following  expression  for  the 
eddy  viscosity: 

^  :z  fs  h  “■  ^ 

k  is  an  empirical  constant,  u  is  the  maximum  velocity  in 
the  mixing  layer  and  is  the  minimum  velocity  in  the 

layer.  The  eddy  viscosity  is  assumed  to  be  constant  across 
the  mixing  layer,  or  a  function  of  only  the  axial  coordinate 
This  Prandtl  eddy  viscosity  model  has  been  widely  applied  to 
a  variety  of  free  mixing  problems  because  of  its  mathemati¬ 
cal  simplicity,  and  because  the  results  agree  fairly  well 
with  experiment  for  many  flow  configurations. 

Methodsof  solution  —  The  two  most  widely  used 
analytical  techniques  for  low  speed  mixing  have  been  (1) 
finding  similarity  solutions  to  the  boundary  layer  equations 
and  (2)  the  von  Kerman  integral  method. 
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The  early  attempts  to  solve  free  mixing  problems 
were  primarily  limited  to  finding  similarity  solutions. 

Using  the  Prandtl  mixing  length  theory  for  turbulent  shear, 
Tollmein  (Ref.  12)  found  solutions  for  (1)  the  two-dimensional 
mixing  layer  between  a  semi— inf inite  stream  and  a  medium  at 
rest,  (2)  the  two-dimensional  jet  flowing  into  a  medium  at 
rest,  and  (3)  the  axisymmetric  jet  flowing  into  a  medium  at 
rest . 

Goertler  (Ref.  13)  applied  the  Prandtl  eddy 
viscosity  model  to  the  two-dimensional  mixing  layer  problem, 
and  obtained  solutions  which  correlated  well  with  the  experi¬ 
mental  results  of  Reichardt  (Ref.  14). 

Early  experiments  on  fully  developed  low  speed 

mixing  zones  and  jets  showed  clearly  that  the  mean  velocity 

profiles  exhibit  shape  similarity,  i.e.,  all  profiles  can 

be  made  approximately  congruent  if  suitably  normalized.  These 

early  experiments  also  indicated  that  the  width  of  the  fully 

developed  mixing  layer  increases  linearly  with  axial  distance, 

X.  This  linear  relation  between  b  and  x  is  inherent  in  the 

analyses  of  Tollmein  and  Goertler.  Tollmein  assumed  that 

i  =  cx,  with  the  proportionality  constant  relating  b  and  x 

being  absorbed  into  the  empirical  constant,  c.  Similarly, 

Goertler  assumed  that  e  =  kx  (u  — u  .  ) . 

'  max  min' 

The  experimental  observation  that  the  fully  de¬ 
veloped  velocity  profiles  are  approximately  similar  has  led 
to  the  widespread  use  of  the  von  Karman  integral  method  for 
the  solution  of  free  mixing  problems.  In  the  integral 
method  the  velocity  profile  shape  is  prescribed,  usually 
based  on  prior  experimental  information,  and  the  boundary 
layer  equations  are  satisfied  only  on  the  average  across 
the  mixing  layer.  It  is  a  well  known  characteristic  of  the 
integral  method  that  the  overall  results  are  relatively  in— 
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sensitive  to  the  profile  shape  chosen,  as  long  as  the  shape 
is  approximately  correct. 

Kuethe  (Ref.  15)  was  the  first  to  use  the  integral 
method  for  the  solution  of  free  mixing  problems.  He  ob¬ 
tained  solutions  for  the  two-dimensional  mixing  layer 
between  two  streams  flowing  at  different  velocities.  Later, 
Squire  and  Trouncer  (Ref.  16)  applied  the  integral  method  to 
the  solution  of  the  axlsymraetrlc  two  stream  mixing  problem. 


The  solutions  of  Kuethe  and  of  Squire  and  Trouncer 
were  obtained  by  using  the  Prandtl  mixing  length  theory  to 
compute  the  turbulent  shear  in  the  mixing  layer.  As  the 
integral  method  is  usually  applied,  with  the  shear  evaluated 
at  one  fixed  point  in  the  normalized  velocity  profile, 

Prandtl 's  mixing  length  and  eddy  viscosity  theories  give 
identical  results.  Both  theories  give  the  same  functional 
form  for  the  characteristic  shear,  i.e.,  proportional  to 

Mixing  with  finite  secondary  velocities  —  It  has 
been  seen  that  Prandtl 's  theories  for  turbulent  shear  have 
been  applied  to  mixing  with  finite  secondary  velocities, 
even  though  the  experiments  which  guided  Prandtl 's  formula¬ 
tion  of  these  theories  were  made  only  with  zero  secondary 
velocity.  Both  theories  predict  that  the  mixing  rate  de¬ 
creases  rapidly  as  the  velocity  ratio, 

increased,  and  that  no  mixing  at  all  occurs  when  =  “niln' 

Experiments,  however,  have  shown  that  the  mixing  rate  of 

species  is  quite  large  when  u  =  u  (Ref.  17).  Fejer, 

mxn 

et  al.  (Ref.  17)  give  an  excellent  discussion  of  this  pro¬ 
blem  of  coaxial  stream  mixing  when  the  secondary  velocity 
is  relatively  large  compared  to  the  primary  velocity. 


Forstall  and  Shapiro  (Ref.  18)  made  a  series  of 
experiments  on  coaxial  stream  mixing  in  the  range  of  secondary- 
primary  velocity  ratios  from  0.2  to  0,5,  It  was  found  that 
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the  analysis  of  Squire  and  Trouncer,  based  on  nixing  length 
theory,  predicted  the  experimental  results  fairly  well  at 
the  lower  velocity  ratios,  but  deviated  considerably  from 
the  experimental  results  at  a  velocity  ratio  of  0.5. 

In  general,  the  Prandtl  theories  for  turbulent 
transport  should  be  applied  with  caution  when  the  secondary- 
primary  velocity  ratio  exceeds  about  0.3. 

Heat  and  mass  transfer  —  Forstall  and  Shapiro 
investigated  the  mass  transfer  In  their  experiments  by  use 
of  a  tracer  gas  In  the  primary  field.  They  also  summarized 
the  results  of  other  investigators  on  heat  and  mass  transfer 
in  free  turbulent  mixing.  In  general,  heat  and  mass  diffuse 
at  the  same  rate  (unity  turbulent  Lewis  number) ,  but  faster 
than  does  momentum.  The  average  value  of  the  turbulent 
Prandtl  (and  Schmidt)  number  is  approximately  0.7  for  low 
speed  mixing. 

2.  Variable  Density  Mixing 

The  Prandtl  theories  for  turbulent  shear  give 
reasonably  good  results  for  nearly  constant  density  mixing, 
as  long  as  the  secondary  velocity  is  small  relative  to  the 
primary  velocity.  When  the  mixing  is  accompanied  by  large 
density  gradients  in  the  layer,  the  current  state  of  know¬ 
ledge  is  much  less  satisfactory. 

Large  transverse  density  gradients  in  the  mixing 
zone  may  be  caused  by: 

1.  Large  difference  in  the  inviscid  stream  Mach 
numbers  (compressibility  effect). 

2.  Large  difference  in  the  temperature  and/or 
compositions  of  the  inviscid  streams. 

3.  Exothermic  chemical  reactions  in  the  mixing 
zones. 
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All  of  these  factors  occur  simultaneously  in  the  most  com¬ 
plex  cases  of  interest. 

The  influence  of  density  gradients  on  the  free 
turbulent  eddy  viscosity  is  currently  a  source  of  speculation. 
Those  who  apply  Prandtl’s  low  speed  shear  theories  directly 
to  variable  density  mixing  problems  assume  that  the  density 
gradient  has  no  influence  on  the  eddy  viscosity  distribution 
in  the  mixing  layer.  At  the  other  extreme  are  those  who  use 
the  compressibility  transformation  technique,  which  implies 
a  large  Influence  of  density  on  the  eddy  viscosity.  The 
reason  for  the  uncertainty  is  that  the  computed  velocity 
profiles  are  relatively  insensitive  to  the  distribution  of 
eddy  viscosity  across  the  mixing  layer.  Consequently,  it  is 
very  difficult  to  determine  the  eddy  viscosity  distributions 
from  mean  flow  measurements  in  variable  density  mixing  ex¬ 
periments.  Unfortunately,  very  precise  measurements  in  a 
mixing  system  with  very  large  density  gradients  have  not  yet 
appeared  in  the  literature. 

Because  of  the  insensitivity  of  the  mean  flow  pro¬ 
files  to  the  eddy  viscosity  model  chosen,  greatly  different 
theories  have  predicted  mean  flow  profiles  which  agree 
fairly  well  with  the  available  experimental  profiles.  The 
empirical  constants  in  the  eddy  viscosity  model,  however, 
usually  must  be  adjusted  to  correlate  the  experimental 
growth  rate  of  the  mixing  layer. 

Methods  of  solution  —  The  various  techniques  for 
analyzing  variable  density  mixing  include  (1)  numerical  so* 
lutions  of  the  boundary  layer  equations,  (2)  compressibility 
transformations,  (3)  linearized  solutions  of  the  boundary 
layer  equations,  and  (4)  integral  methods. 

Vasiliu  (Ref.  19)  obtained  numerical  solutions  for 
the  problem  of  the  mixing  and  combustion  layer  along  the 
boundary  between  a  supersonic  rocket  plume  and  surrounding 
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air  stream.  The  axial  pressure  gradient  was  prescribed  from 
a  prior  inviscid  analysis,  and  finite  rate  chemistry  was  in¬ 
cluded  in  the  mixing  analysis.  The  eddy  viscosity  was  assumed 
to  be  a  function  of  only  the  streamwise  coordinate,  and  non- 
unity  (but  constant)  turbulent  Prandtl  and  Lewis  numbers 
could  be  used.  Such  solutions  are  a  tribute  to  modern  cort- 
puter  technology,  but  unfortunately  the  current  uncertainty 
about  turbulent  transport  properties  makes  it  doubtful 
whether  this  approach  is  any  more  rewarding  than  simpler 
methods.  The  current  knowledge  about  chemical  reaction  rates 
in  turbulent  flows  is  even  poorer  than  that  about  turbulent 
transport.  After  the  knowledge  about  turbulent  transport 
and  chemical  reaction  rates  has  progressed,  the  numerical 
technique  will  undoubtedly  become  the  most  powerful  approach 
for  free  mixing  problems. 

Mager  (Ref.  20)  proposed  a  compressibility  trans¬ 
formation  of  the  boundary  layer  equations  which  is  the 
turbulent  counterpart  of  the  Howarth  transformation  for 
laminar  flow.  Ting  and  Libby  (Ref.  21)  applied  Mager 's 
transformation  to  free  turbulent  mixing  and  derived  an  ex¬ 
pression  for  the  eddy  viscosity.  Libby  (Pef.  22)  analyzed 
isobaric  coaxial  jet  mixing,  including  equilibrium 
chemical  reactions,  with  the  ccmipressibility  transformation 
technique.  The  corresponding  incompressible  solutions  were 
obtained  from  linearized  equations  incorporating  the  Prandtl 
eddy  viscosity  model.  The  turbulent  Prandtl  and  Lewis  num¬ 
bers  were  assumed  to  be  unity.  A  reference  density  appears 
in  the  transformation,  and  this  must  be  chosen  arbitrarily 
because  the  theory  offers  no  guidance  on  this  subject. 

The  compressibility  transformation  predicts  the 
largest  distortion  of  the  incompressible  velocity  profile 
shape  of  any  of  the  methods  in  current  use.  The  available 
experimental  results,  although  limited,  do  not  tend  to  sub— 
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stantiate  the  profile  distortion  predicted  by  the  compressi¬ 
bility  transformation  (Ref.  23,  for  example).  In  fact, 
Peters  (Ref.  24)  used  this  experimental  observation  that 
the  velocity  profile  shape  is  apparently  invariant  to 
develop  a  speculative  model  for  the  variable  density  eddy 
viscosity. 

Kleinstein  (Ref.  25)  presented  a  linearized  solu¬ 
tion  of  the  compressible  boundary  layer  equations  for 
isobar ic  coaxial  stream  mixing.  He  assumed  that  the  product 
pG  Is  a  function  of  only  the  axial  coordinate.  One  is  in¬ 
clined  to  prefer  Kleinstein 's  solution  over  the  compressi¬ 
bility  transformation  solution  because  it  predicts  less 
distortion  of  the  incompressible  profile  shapes. 

Because  the  mixing  zone  velocity  profiles  are  in¬ 
sensitive  to  transverse  density  gradients,  the  integral 
method  is  useful  for  analyzing  variable  density  mixing 
problems.  The  normalized  profile  shape  is  usually  assumed 
to  be  the  same  as  the  constant  density  profile  shape.  One 
advantage  of  the  integral  method  is  that  it  is  relatively 
easy  to  include  the  effects  of  streamwise  pressure  gradients 

Korst,  et  al.  (Ref.  26)  have  used  the  integral 
technique  quite  successfully  to  analyze  the  two-dimensional 
compressible  mixing  layer,  with  application  to  the  flow  past 
back  steps.  Warren  (Ref.  27)  and  Donaldson  and  Gray  (Ref. 
28)  have  used  the  integral  method  to  analyze  the  variable 
density  mixing  of  an  axisymmetric  Jet  with  a  quiescent 
medium.  Donaldson  and  Gray  also  proposed  a  simple  cor»- 
pressibility  correction  to  the  Prandtl  eddy  viscosity  model, 
based  on  a  characteristic  Mach  number  in  the  mixing  layer. 

Heat  and  mass  transfer  —  As  is  the  case  for  the 
eddy  viscosity,  the  relative  rates  of  transport  for  heat, 
mass  and  momentum  are  not  well  understood  for  variable 


19 


AEDC-TR-68.270 


density  mixing.  Alpinieri  (Ref.  29)  and  Zakkay,  et  al. 

(Ref.  30)  made  a  series  of  experiments  on  nonreactive  co¬ 
axial  stream  mixing.  The  outer  stream  was  either  subsonic 
or  supersonic  air,  and  the  jet  fluid  was  either  hydrogen, 
helium,  argon  or  carbon  dioxide.  In  addition  to  showing  that 
the  mixing  rate  is  not  zero  when  the  streams  have  equal 
velocities,  these  experiments  Indicated  that  heat  and  mass 
tend  to  diffuse  faster  than  momentum. 

Zakkay,  et  al.  evaluated  the  centerline  values  of 
the  turbulent  Schmidt  and  Lewis  numbers  for  their  various 
experiments,  and  found  a  great  deal  of  scatter  in  the  re¬ 
sults.  Most  of  the  Schmidt  numbers  fell  between  0.4  and  2.0. 
Somewhat  less  scatter  was  observed  for  the  Lewis  number 
results,  with  most  of  the  data  falling  between  0.5  and  1.2. 
These  results  show  that  it  is  extremely  questionable  whether 
the  relative  transport  of  momentum,  heat,  and  mass  can  be 
accurately  represented  by  constant  values  of  the  Prandtl 
(or  Schmidt)  and  Lewis  numbers.  No  experimental  results 
have  come  to  the  attention  of  the  author  on  the  relative 
rates  of  transport  of  momentum,  heat  and  mass  for  mixing 
with  simultaneous  chemical  reactions. 

3.  Ducteii  Coaxial  Stream  Mixing 

When  the  two  streams  mix  Inside  a  duct,  the  axial 
pressure  gradient  is  not  prescribed,  but  must  be  computed 
as  an  additional  unknown  parameter  in  the  flow  field  solu¬ 
tion. 

Mikhail  (Ref.  31)  used  the  integral  method  to 
analyze  the  low  speed  flow  in  a  constant  area  jet  pump.  He 
assumed  that  the  mixing  zone  velocity  profile  could  be 
represented  by  a  cosine  function,  and  used  Prandtl 's  mixing 
length  theory  to  evaluate  the  shear  along  one  control  sur¬ 
face  in  the  mixing  zone. 
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For  ducted  mixing  with  large  density  gradients,  re¬ 
latively  few  analyses  have  come  to  the  attention  of  the  author.* 
Emmons  (Ref.  32)  treated  the  flow  in  the  mixing  duct  of  an  air 
augmented  rocket  by  numerical  solution  of  the  boundary  layer 
equations.  Finite  rate  chemistry  was  considered  in  the  mixing 
layer,  although  with  a  highly  simplified  reaction  model.  The 
eddy  viscosity  was  assumed  to  be  a  function  of  only  the  stream- 
wise  coordinate.  Emmons'  paper  is  very  brief  and  it  is  not 
clear  how  he  treated  the  axial  pressure  gradients.  He  ap¬ 
parently  assumed  the  inviscid  portions  of  the  flow  to  be  one- 
dimensional,  although  this  point  is  also  unclear. 

Peters,  et  al.  (Ref.  6)  used  the  integral  method 
to  analyze  ducted  coaxial  stream  mixing  with  equilibrium 
chemical  reactions.  The  turbulent  Prandtl  and  Lewis  numbers 
were  assumed  to  be  unity,  and  the  normalized  mixing  zone 
velocity  profile  was  assumed  to  be  a  cosine  function.  Rather 
than  characterizing  the  growth  rate  of  the  mixing  zone  with 
a  shear  stress  evaluated  along  a  control  surface  in  the  zone, 
the  equation  of  Abramovich  (Ref.  33)  was  used  to  calculate 
the  growth  rate.  The  assumption  that  the  inviscid  portions 
of  the  flow  were  one— dimensional  limited  the  validity  and 
applicability  of  the  theory  for  supersonic  mixing,  although 
the  predicted  mixing  zone  profiles  agreed  fairly  well  with 
experimental  profiles.  Another  limitation  of  the  theory  was 
the  use  of  the  Abramovich  mixing  zone  growth  rate  equation, 
which  is  strictly  applicable  only  to  isobar ic  mixing. 

The  analysis  presented  in  this  report  may  be  con¬ 
sidered  an  extension  and  refinement  of  the  theory  of  Ref.  6. 


*  Edelman  and  Fortune  (Ref.  47)  recently  presented  a  boundary 
layer  type  analysis  (negligible  lateral  pressure  gradients) 
for  ducted  coaxial  stream  mixing.  The  viscous  layer  is  com¬ 
puted  by  finite  difference  methods,  and  finite  rate  chemistry 
is  included  in  the  computations.  Either  the  duct  pressure 
distribution  or  the  duct  shape  can  be  prescribed. 
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4.  Ducted  Mixing  With  Nenuniform  Inviscid  Core 

For  certain  classes  of  problems  in  ejectors  and 
propulsion  systems,  it  is  essential  that  the  nonuniform 
nature  of  the  inviscid  core  flow- be  considered.  For  example, 
a  supersonic  ejector  system  having  a  conical  nozzle  and 
operating  in  the  "upstream  choking"  mode  (Section  I)  cannot 
be  satisfactorily  analyzed  if  the  inviscid  core  flow  is 
assumed  to  be  one— dimensional .  The  amount  of  pluming  of  the 
supersonic  primary  stream  depends  on  the  detailed  flow  pro¬ 
file  at  the  primary  nozzle  exit. 

Addy  (Ref.  34)  and  Chow  and  Addy  (Ref.  5)  have 
analyzed  the  flow  in  an  axisynimetric  ejector  system.  The 
secondary  flow  was  assumed  to  be  one— dimensional  and  the 
primary  stream  was  treated  with  the  method  of  characteristics. 
The  inviscid  flow  field  was  first  computed,  then  the  effects 
of  mixing  were  obtained  by  superimposing  the  two-dimensional 
mixing  profiles  of  Korst  and  Chow  (Ref.  35)  on  the  inviscid 
boundary.  Chow  and  Addy  obtained  excellent  agreement  between 
their  analysis  and  their  air— air  ejector  experiments.  The 
superposition  technique  gives  good  results  for  the  air— air 
ejector  operating  in  the  upstream  choking  mode  because  the 
effect  of  mixing  is  relatively  small*. 

Chow  and  Addy  clearly  point  out  that  their  solu¬ 
tion  is  valid  only  for  flows  with  relatively  weak  viscous 
interaction.  They  suggest  that  their  solution  can  be  im¬ 
proved  by  successive  viscous  corrections  to  the  inviscid 


*  The  discussion  here  is  limited  to  cases  in  which  the  secon¬ 
dary  flow  is  large  enough  so  that  the  mixing  zone  does  not 
impinge  on  the  duct  wall.  For  zero  and  low  flow  rates,  the 
mixing  zone  impinges  on  the  wall  and  viscous  effects  are 
very  important.  Chow  and  Addy  analyzed  these  low  secondary 
flow  cases  by  using  a  base  pressure  analysis  which  included 
effects  of  base  bleed. 
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solution,  but  that  the  superposition  technique  is  limited  to 
flows  with  thin  mixing  layers,  i.e.,  approximately  two- 
dimensional  mixing.  Another  limitation  of  the  superposition 
technique  in  general  is  that  the  mixing  profiles  are  com¬ 
puted  with  the  assumption  that  there  are  no  transverse 
pressure  gradients,  but  then  are  superimposed  onto  a  flow 
which  has  strong  radial  pressure  gradients. 

To  be  valid  in  the  downstream  portion  of  the  flow 
field,  as  well  as  to  be  self  consistent  in  the  upstream 
portion,  the  inviscid  flow  field  and  the  mixing  layer  must 
be  computed  simultaneously.  The  simultaneous  solution  of 
the  inviscid  core  flow  (with  the  method  of  characteristics) , 
the  mixing  layer  and  the  secondary  inviscid  flow  is  the 
subject  of  Part  II  of  this  report  (2— D  Core  Theory) . 
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SECTION  III 

GENERAL  DISCUSSION  OF  METHODS  USED  IN  THE  ANALYSIS 

The  purpose  of  this  section  is  to  outline  the 
principal  assumptions  and  techniques  used  in  this  investiga¬ 
tion. 

1.  Outline  of  Principal  Assumptions 

As  pointed  out  in  Section  II,  the  analysis  of  free 
turbulent  mixing  is  necessarily  approximate,  and  a 
number  of  assumptions  is  required  to  provide  a  complete 
mathematical  representation  of  the  problem.  These  as¬ 
sumptions,  as  well  as  others  which  have  been  made  for 
mathematical  convenience,  are  listed  below.  The  more 
critical  or  questionable  assumptions  are  discussed 
briefly. 

(1)  The  flow  is  axisymmetric, 

(2)  All  gases  obey  the  perfect  gas  law. 

(3)  The  usual  assumptions  used  in  the  analysis  of 
turbulent  boundary  layers  are  assumed  to  be  appli¬ 
cable.  The  most  important  of  these  assumptions  is 
that  the  radial  pressure  gradients  in  the  mixing 
layer  are  negligible.  The  neglect  of  radial  pres¬ 
sure  gradients  has  yielded  good  results  for  low 
speed  free  mixing  problems  (Ref.  8),  hut  may  not 
be  valid  when  the  mixing  layer  is  wholly  or  partly 
supersonic.  Nevertheless,  this  assumpxion  has  been 
widely  applied  in  the  analysis  of  high  speed  mixing 
layers.  It  is  not  intended  that  the  neglect  of 
radial  pressure  gradients  be  justified,  but  rather 
that  the  assumption  be  clearly  noted, 

(4)  The  inviscid  portion  of  the  secondary  flow  is 
assumed  to  be  one— dimensional  and  isentropic.  This 
assumption  is  quite  valid  when  the  secondary  flow 
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Is  subsonic  and  the  radii  of  curvature  of  the 
streamlines  are  large  with  respect  to  the  radius 
of  the  flow  passage.  This  mixing  theory  is  also 
applied  to  cases  in  which  the  secondary  flow 
accelerates  into  the  low  supersonic  range  (Ma  <  1.3). 
The  one— dimensional  assumption  will  be  less  valid 
for  these  supersonic  cases,  but  this  is  not  con¬ 
sidered  a  serious  defect  in  the  analysis  compared 
to  some  of  the  more  questionable  assumptions. 

(5)  The  invlscid  portion  of  the  primary  flow  (Inviscid 
core)  is  assumed  to  be  one— dimensional  and  isen— 
tropic  for  the  1— D  Core  Theory  (Part  I) .  Again, 
this  assumption  will  be  quite  valid  for  subsonic 
primary  flows.  For  supersonic  primary  flows,  this 
assumption  will  be  approximately  correct  only  if 
the  primary  flow  is  initially  uniform  and  the  whole 
process  is  nearly  isobar ic. 

(6)  The  inviscid  core  flow  is  assumed  to  be  isen— 
tropic  for  the  2-D  Core  Theory  (Part  II).  This 
region  of  the  flow  is  treated  with  the  irrota— 
tlonal  method  of  characteristics.  For  primary 
nozzles  which  provide  uniform  or  conical  flow,  it 
is  well  known  that  the  shock  waves  in  the  free  Jet 
are  weak  in  the  near  field,  but  become  more  signi¬ 
ficant  in  the  far  field  (several  primary  nozzle 
radii  downstream).  Some  of  the  mixing  configurations 
considered  in  this  investigation  are  relatively  long, 
and  it  is  to  be  expected  that  neglecting  the  shock 
waves  in  the  flow  will  cause  some  deviation  from 
experiment . 

(7)  The  mixing  layer  is  assumed  to  be  completely  tur¬ 
bulent,  and  the  boundary  layer  thicknesses  at  the 
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initiation  of  the  free  mixing  layer  are  assumed 
to  be  negligible.  The  first  of  these  assumptions 
implies  that  transition  from  laminar  to  turbulent 
flow  occurs  well  upstream  of  the  initial  section 
of  the  mixing  duct.  The  assumption  that  the 
initial  boundary  layers  are  thin  is  reasonably 
valid  in  this  analysis  because  the  major  interest 
is  in  the  region  many  initial  boundary  layer 
thicknesses  downstream  from  the  mixing  duct  entrance. 

(8)  The  finite  thickness  of  the  nozzle  lip  separating 
the  primary  and  secondary  flows  is  assumed  to  be 
negligible. 

(9)  The  viscous  effects  at  the  duct  wall  are  assumed 
to  be  negligible.  In  one  correlation  of  the  flow 
in  a  low  speed  jet  pump  (Section  IX)  the  duct  is 
approximately  10  diameters  long  and  the  wall 
boundary  layer  effects  are  significant.  In  that 
one  case  the  wall  viscous  effects  were  computed 
by  use  of  an  approximate  displacement  thickness 
correction  to  the  duct  contour. 

(10)  The  shape  of  the  normalized  velocity  profiles  in 
the  mixing  layer  is  assumed  to  be  similar  at  all 
axial  locations,  and  to  be  represented  by  a  cosine 
function. 

(11)  The  turbulent  Prandtl  and  Lewis  numbers  in  the 
mixing  layer  are  assumed  to  be  unity. 

(12)  For  mixing  with  simultaneous  chemical  I'eactions, 

the  reactions  are  assumed  to  be  in  equilibrium.  The 
current  knowledge  of  reaction  chemistry  in  turbulent 
flows  is  insufficient  to  allow  a  theoretical  verifi¬ 
cation  of  this  assumption;  its  validity  can  only  be 
checked  by  comparison  of  the  theory  with  experiment. 
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(13)  The  eddy  viscosity  at  the  half-velocity  point  in 
the  mixing  layer  is  assumed  to  be  given  by  the 
Prandtl  eddy  viscosity  equation,  with  the  effect 
of  variable  density  treated  by  a  modification  of 
the  correction  suggested  by  Donaldson  and  Gray. 

This  eddy  viscosity  model  is  discussed  further  in 
Section  VI. 

(14)  In  the  third  regime,  where  the  turbulent  region 
extends  across  the  entire  duct  (Fig,  1),  the  free 
mixing  concepts  of  shear  and  profile  shape  simi¬ 
larity  are  assumed  to  be  applicable.  Mikhail 
(Ref.  31)  found  that  the  cosine  profile  was  a 
reasonably  good  representation  of  the  third  re¬ 
gime  velocity  profiles  in  his  low  speed  Jet  pump 
configuration.  It  is  realized  that  the  third 
regime  flow  will  eventually  become  a  fully  de¬ 
veloped  duct  flow  in  which  wall  effects  pre¬ 
dominate.  However,  this  analysis  is  not  intended 
for  application  to  more  than  1  or  2  duct  diameters 
downstream  into  the  third  regime.  The  validity 

of  free  turbulent  shear  relations  in  the  upstream 
portion  of  the  third  regime  can  only  be  established 
by  comparison  of  theory  and  experiment. 

These  assumptions,  along  with  others  that  are  mentioned 
in  the  text,  sufficiently  define  the  ducted  mixing  pro¬ 
blem  so  that  flow  field  solutions  can  be  obtained  by 
application  of  the  integral  method. 

2.  Outline  of  Analytical  Developments  -  1-D  Core  Theory 

Section  IV  presents  the  derivation  of  three  basic 
integral  equations  for  the  ducted  mixing  problem. 

These  integral  equations,  which  are  derived  from  the 


27 


AEDC-TR-68-270 


continuity  equation  and  the  axial  boundary  layer  momentum 
equation,  are: 

(1)  A  continuity  equation  for  the  entire  flow. 

(2)  A  momentum  equation  for  the  entire  flow. 

(3)  A  momentum  equation  for  the  portion  of  the  flow  be¬ 
tween  the  duct  centerline  and  a  control  surface 
located  halfway  across  the  mixing  layer. 

A  term  representing  the  turbulent  shear  stress  appears 
in  the  last  equation. 

Section  V  presents  the  equations  for  the  mixing  zone 
velocity,  and  for  relating  the  mixing  zone  temperature 
and  density  to  the  velocity.  The  method  of  treating  the 
mixing  zone  chemistry  is  presented. 

Section  VI  presents  the  model  used  for  the  turbulent 
eddy  viscosity.  A  modification  of  the  density  correction 
suggested  by  Donaldson  and  Gray  is  developed.  The  range 
of  applicability  of  the  eddy  viscosity  model  is  discussed. 

Section  VII  presents  the  method  of  transforming  the  set 
of  three  basic  differential— integral  equations  (derived 
in  Section  IV)  into  a  set  of  three  first  order  ordinary 
differential  equations  for  three  variables  which  fully 
define  the  flow.  All  of  the  auxiliary  functions  which 
are  needed  to  compute  the  coefficients  of  these  dif¬ 
ferential  equations  are  presented  in  Appendix  I.  The 
resulting  system  of  equations  is  solved  by  standard 
numerical  techniques. 

Section  VIII  presents  a  brief  discussion  of  the  numerical 
aspects  of  the  solution,  from  the  point  of  view  of  the 
engineering  results. 
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3.  Outlin«  of  Analytical  Developments  •>  2*D  Core  Theory 

Section  X  presents  the  derivation  of  three  basic  Integral 
equations.  These  equations  differ  from  the  1— D  Core 
Theory  equations  in  that  the  Integrals  do  not  include 
the  inviscid  core;  the  influence  of  the  inviscid  core 
appears  in  terms  of  the  boundary  conditions  for  the 
integral  equations. 

Section  XI  presents  the  equations  used  to  represent  the 
flow  in  the  mixing  layer.  With  a  few  exceptions,  the 
mixing  zone  parameters  are  determined  in  the  same  way 
as  in  the  1— D  Core  Theory. 

Section  XII  presents  the  transformation  of  the  three 
basic  differential-integral  equations  (derived  in  Section 
X)  into  a  set  of  ordinary  differential  equations  suitable 
for  numerical  solution  by  conventional  techniques.  All 
of  the  auxiliary  functions  required  to  compute  the 
coefficients  are  presented  in  Appendix  H. 

Section  XIII  presents  the  method  of  characteristics 
solution  for  the  inviscid  core  flow.  The  key  technique 
for  coupling  the  method  of  characteristics  solution  to 
the  system  of  integral  equations  is  discussed. 

Section  XIV  presents  a  brief  discussion  of  the  numerical 
aspects  of  the  2— D  Core  Theory  solutions. 

Section  XV  presents  the  methods  used  to  obtain  solutions 
for  ducted  mixing  systems  operating  in  (1)  the  "down¬ 
stream  choking"  mode,  (2)  the  "back  pressure  dependent" 
mode,  and  (3)  the  "upstream  choking"  mode.  In  all  cases, 
the  solutions  are  obtained  by  iteration  of  the  initial 
static  pressure. 

Appendix  III  presents  a  listing  of  computer  program  in¬ 
puts  and  printout  variables  for  both  the  1— D  Core  Theory 
and  2— D  Core  Theory  programs. 
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PART  I 

DUCTED  MIXING  WITH  ONE-DIMENSIONAL  INVISCID  CORE  FLOW 
(1-D  CORE  THEORY) 

SECTION  IV 

DEVELOPMENT  OF  BASIC  DIFFERENTIAL  -  INTEGRAL  EQUATIONS 

In  this  section  a  set  of  three  basic  differential— 
integral  equations  is  derived  from  the  boundary  layer  equa¬ 
tions.  The  three  equations  are  (1)  a  continuity  equation  for 
the  entire  flow,  (2)  a  momentum  equation  for  the  entire  flow, 
and  (3)  a  momentum  equation  for  the  portion  of  the  flow  between 
the  duct  centerline  and  an  arbitrary  control  surface. 

1.  Integration  of  Boundary  Layer  Equations 

The  axial  momentum  equation  for  axisyminetr ic  flow  is: 

H  ^  It  =  H;  (1) 

The  continuity  equation  is: 

Ij,  Cru.f)  4-  =  o  (2) 

The  pressure  is  a  function  of  only  the  axial  coordinate, 

and  is  equal  to  the  wall  pressure,  p  .  After  multiplying 

w 

each  term  by  dr,  eq.  (1)  is  integrated  between  the  ar¬ 
bitrarily  selected  lower  and  upper  limits,  r^  and  r2: 


The  term  "dif ferential— integral  equations"  is  used  to 
denote  the  integral  form  of  the  boundary  layer  differential 
equations.  These  equations  contain  terms  of  the  type 


1  r  ^ 
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Now  consider  the  second  term  of  eq.  (3).  Integrating  by 
parts,  ^  fa 

f?vT~  Ar  =  -  u,  |-(fi;-'r)  At 

'T  ' 

'  1  ^  I  'T 

^  I 

Noting  that  from  eq.  (2),  eq.  (4) 

is  substituted  into  eq.  (3):  ^ 


(4) 


u  eii.  V.J7  ; 

PtU^a^s-^  ~  f»,Lr;  ai'f;  +  Ju  d'T  =i 


dxi  j 


Tk. 


Combining  terms: 

nr.-r,^-  ^CV 


or 

'T* 


V  civ  -H  Px.^2,>-^2.'<”z  -  Uaj  'Cl  =. 

-  ^ 
d^pC 


Vci'C 


(5) 


I 

Now  let  r-  =  0  and  r^  =  r  ,  and  note  that  t„  =  0: 

X  u  w  w 

Ty.  ^ 

[^(pa^)TA'<-  -(-  P^LruUy^'Cyi  =.  -  ^ 

<iK  z- 


(6> 
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Integrating  the  continuity  equation  (eq.  2) 
pu-'T  =  - 
o 


Thus 


a'T 

o 

Substituting  eq.  (7)  into  eq.  (6): 

O  ' 

O 

To  provide  a  second  momentum  equation,  eq.  (5)  is 
evaluated  with  =  0  and  V2  =  the  arbitrary  position 


in  the  mixing  zone  where  u  =  u  : 

HI 


JOo  .'Qr> 

-a„  j|- (ftf,>a>r=  ^ 

3  O  ^ 


(7) 


(8) 


(9) 


Equations  (7) ,  (8) ,  and  (9)  are  the  three  basic 

integral  equations,  a  continuity  equation  for  the 

entire  flow,  a  momentum  equation  for  the  entire  flow  and 

a  momentum  equation  for  the  interior  portion  of  the  flow 

between  r  =  0  and  r  =  r  . 

m 

Applying  Leibnitz’  rule  for  differentiating  an 
integral  (Ref.  36)  to  eq.  (7): 

O 


rw 


fit'Tci'r  ■*“  Pw  '^Vyi 


Because  the  wall  is  a  streamline  in  the  flow  field, 

— ^  continuity  equation  becomes; 
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d 


f 


?OL'rAT 


Q 


-  o 


Similarly,  applying  Leibnitz'  rule  to  eq.  (8)  gives: 


or 


d  t 


rcir  = 


_  d^  r*,* 


d;d 


Applying  Leibnitz'  rule  to  eq.  (9)  gives 

wi 


'^rrk 


d  Z 


_a,^^jp.64d«r 


4*  Uinfl  'f'm  ■4"^')  ~  “*  ^^2. 

O^c-  /  d  )4  ^ 


^  '  TTI 

i  j PM.*-r^r - 

o  o 

Equations  (10),  (11)  and  (12)  are  the  working  forms  of  the 
three  basic  integral  equations. 


(10) 


(11) 


(12) 
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2.  First  Regime  Differential-Integral  Equations 

Referring  to  the  schematic  shown  in  Fig.  3,  the 
mixing  layer  in  the  first  regime  is  bounded  on  either  side 
by  a  region  of  one—dimensional  isentropic  flow.  The  radius 
of  the  mixing  zone  control  surface,  r^^,  is  taken  to  be  half 
the  radial  distance  across  the  mixing  zone,  i.e.,  r^  = 
r^  +  b/2.  The  integrals  appearing  in  eqs.  (10),  (11)  and 
(12)  can  be  divided  into  parts  corresponding  to  the  various 
regions  of  the  flow.  Thus  the  following  equations  are 
obtained. 


Continuity  Equation 

,'C:-hh 


" 


(13) 


Total  Momentum  Equation 


_  _  dfw  SV 
z 


(14) 


Half— Radius  Momentum  Equation 

S;  + 


2. 

-  pardr  =  r„('ri+  t/a.)  - 

^  d%  Z 


(15) 
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3.  Second  Regime  Di'fferential-lntegrol  Equations 

In  the  second  regime  (Fig.  3)  the  mixing  zone 
extends  to  the  centerline  (r^  =  0) ,  but  a  region  of  one- 
dimensional  isentropic  flow  exists  outside  the  mixing 
layer.  In  this  regime,  r^^  =  b/2,  and  the  following  equa<- 
tions  are  obtained. 


Continuity  Equation 

■h 


d  r 


d;t 


=  O 


Total  Momentum  Equation 

o 

Half— Radius  Momentum  Equation 


- 

A7L  Z 


i  (pu.Vd'r-Um^ifitTcl'r  =.tV^4  -  dfi  ^ 
dx.)  ^  djA  s 


4.  Third  Regime  Differential— Integral  Equations 


In  the  third  regime  (Fig-  3)  the  mixing  zone 
extends  across  the  entire  duct,  thus  r.  =0,  b  =  r  and 
~  ^w/2*  integral  terms  in  eqs.  (10),  (11)  and  (12) 

only  include  the  mixing  layer,  and  the  following  equations 
are  obtained. 


(16) 


(17) 


(18) 


Continuity  Equation 


A. 


fa'T 


o 


o 


(19) 


36 


AEDC‘TR-68-270 


Total  Momentum  Equation 


ax- 


o 


Half—Radlus  Momentum  Equation 


(20) 


(21) 
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SECTION  V 

MIXING  ZONE  PROFILES 


A  set  of  three  basic  differential— integral  equa¬ 
tions  was  developed  in  Section  IV  for  each  of  the  three 
mixing  regimes.  Solution  of  these  differential-integral 
equations  requires  that  p  and  u  in  the  nixing  layer  be  de¬ 
fined  in  terms  of  the  other  variables  of  the  problem.  The 
following  development  of  the  relations  for  the  density  and 
velocity  distributions  is  essentially  the  same  as  that  of 
Peters,  et  al.  (Ref.  6). 

1.  Relation  Between  Velocity,  Composition  end  Enthalpy 

The  axisymmetric  boundary  layer  equations  for 
momentum,  energy  and  conservation  of  elemental  species,  with 
unity  turbulent  Prandtl  and  Lewis  numbers,  are  (Ref.  22): 

Momentum  Equation 


Energy  Equation 

p«.r  iib  + 


Ho 


Conservation  of  Elemental  Species 


f  ar 


pir>r 


These  equations  are  the  steady  laminar  boundary  layer 
equations  with  the  laminar  transport  properties  replaced  by 
the  corresponding  effective  turbulent  values.  For  free 
turbulent  mixing,  the  turbulent  transport  is  normally  an 
order  of  magnitude  or  moi'e  larger  than  the  molecular 


(22) 


(23) 


(24) 
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transport;  therefore,  the  molecular  transport  coefficients 
are  neglected. 

For  application  to  flows  with  chemical  reaction, 

It  is  Important  to  note  that  is  defined  as  the  local 

mass  fraction  of  the  element  k,  and  is  defined  as  the 

stagnation  enthalpy,  including  chemical  heats  of  formation. 

As  Libby  has  shown  (Ref.  22),  these  definitions  for  and 

H  remove  the  usual  formation  terms  from  the  energy  and 
o 

species  equations.  For  chemically  frozen  flow,  can  be 
Interpreted  as  the  local  mass  fraction  of  the  molecular 
species  k,  and  can  be  defined  with  the  usual  definition 
of  enthalpy,  l.e., 

o  ^ 

Now  consider  the  case  of  isobaric  mixing  of  an 
initially  uniform  jet  with  a  surrounding  uniform  infinite 
stream.  The  dp/dx  term  In  eq.  (22)  is  zero,  and  eqs.  (22), 
(23)  and  (24)  are  identical  in  form.  The  Initial  and 
boundary  conditions  for  all  variables  are  similar,  and 
the  following  linear  relation  is  obtained  between  u, 
and  Cj^: 

The  subscript  ’’a”  refers  to  the  outer  stream  and  "j"  to 
the  inviscid  jet.  and  Cj^  are  related  to  u  by  eq.  (25) 
for  isobaric  mixing,  and  if  these  variables  also  define 
the  density,  p,  then  the  solution  of  the  momentum  equation 
completely  defines  the  flow  field. 

Equation  (25)  is  also  assumed  to  be  valid  for 
mixing  with  axial  pressure  gradients.  Note,  however,  that 
the  Inviscid  reference  velocities,  u^  and  Uj ,  are  pressure 
dependent.  The  assumption  that  eq.  (25)  is  valid  for  flows 
with  axial  pressure  gradients  is  made  conditionally.  If 


(25) 


39 


AEDC-TR-68-270 


eq.  (25)  is  valid,  then  the  computed  total  flux  of  any 
elemental  species  should  be  the  same  at  all  axial  stations  in 
the  flow  field.  In  the  ducted  mixing  case,  this  species 
conservation  may  be  expressed  as 


ns 


■tcur\'h 


If  the  integral  of  eq.  (26)  does  remain  approximately 
constant  throughout  a  flow  field  computation,  then  the 
use  of  eq.  (25)  is  justified,  at  least  for  an  Integral 
method.  In  other  words,  the  solution  of  the  integrated 
momentum  equation  is  also  a  satisfactory  solution  of  the 
integrated  species  and  energy  equations. 


(26) 


Because  the  transport  coefficients  for  all 
species  are  the  same,  the  composition  at  any  point  in  the 
mixing  zone  may  be  considered  a  simple  two  component  mix¬ 
ture.  This  composition  may  be  characterized  by  the  single 
parameter  which  is  defined  as  the  mass  fraction  of 
elements  originating  in  the  central  stream  (thus  U.  =  1, 

C  =  0) .  Equation  (25)  becomes 

A 


CL  “•  U0c^  /Vo  —  Hocu  ^  ^ 

Considering  eq.  (26),  a  species  conservation 
parameter,  Q,  is  now  defined  in  terms  of  C, 


where  w.  is  the  initial  mass  flow  of  the  jet.  If  eq.  (25a) 

XJ 

is  correct,  then  the  computed  value  of  Q  will  be  equal  to 
unity  throughout  the  flow  field. 


(25a) 


(27) 


The  theoretical  variation  of  Q  for  a  ducted  rocket- 
air  mixing  system  is  discussed  in  Section  XVII. 
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2-  Mixing  Zon»  Velocity  Profiles 


It  is  assumed  that  the  mixing  zone  velocity  pro¬ 
files  are  similar  and  can  be  represented  by  a  cosine  function: 


CC  ml  fk  ^ 


(28) 


In  the  first  regime  (Fig.  3)  the  velocity  profile  equation 
is 


a  - 

LL^-Ue^ 


(28a) 


It  has  been  assumed  that  the  radial  pressure  gradients  In 

the  mixing  zone  are  negligible,  and  that  the  inviscld 

secondary  and  primary  flows  are  one-dimensional.  It  follows 

then  that  the  pressure  across  the  entire  duct  at  any  axial 

station  is  p^,  with  p^  being  a  function  only  of  x.  The 

inviscld  stream  velocities,  Uj  and  u^,  are  isentropic 

functions  of  only  p^,  therefore  eq.  (28a)  defines  u  In  terms 

of  p  ,  r . ,  and  b. 
w  1 


is 


In  the  second  regime  the  velocity  profile  equation 


CCou 


** 

'4  I  +  cos  TrC-f-) 

L 


(28b) 


where  u  is  the  centerline  velocity  (Fig.  3) .  In  the  third 
c 

regime  the  velocity  profile  equation  is 

=  4  [f  ^ 

L-  J 


(28c) 


where  u„  is  the  velocity  at  the  duct  wall, 
w 
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Note  that  the  composition  profiles  are  determined 
from  these  velocity  profile  equations  (eqs.  28a-c)  and  eq. 
(25a) . 

Velocity  along  the  half— radius  control  surface  — 

The  cosine  velocity  profile  gives  the  result  that  the  velocity 
at  the  half-radius  point  in  the  profile  is  the  average  of 
the  inviscid  stream  velocities.  In  the  first  regime  the 
velocity,  u^^^,  is  given  by 

In  the  second  regime 

In  the  third  regime 


(29a) 

(29b) 

(29c) 


3.  Mixing  Zone  Deniity 

The  mixing  25one  density  is  given  by  the  perfect 

gas  law, 

P  ■=  ^/RT  (30) 

where  T  is  the  static  temperature  and  R  is  the  gas  constant, 
defined  as  the  universal  gas  constant  divided  by  the  molecular 
weight  of  the  gas. 

The  method  cf  determining  the  mixing  zone  tempera¬ 
tures  and  gas  constants  is  discussed  in  the  next  two  sections 
for  the  cases  of  equilibrium  and  frozen  chemistry,  respec¬ 
tively. 

4.  Equilibrium  Mixing  Zone  Chemistry 

For  rocket— air  mixing,  the  mixing  zone  chemistry  is 
assumed  to  be  in  equilibrium,  i.e.,  the  chemical  reaction 
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rates  are  assumed  to  be  very  fast  with  respect  to  a  charac¬ 
teristic  time  for  the  flow.  The  temperatures  and  gas 
properties  in  the  mixing  zone  must  be  computed  by  application 
of  the  concepts  of  equilibrium  chemistry.  If  the  mixing  were 
Isobaric,  then  the  mixing  zone  temperatures  and  gas  proper¬ 
ties  would  be  unique  functions  of  T!  for  any  given  problem. 

In  the  case  of  ducted  mixing,  the  pressure  will  not  be  con¬ 
stant  in  any  given  problem,  and  in  principle,  the  mixing 
zone  temperatures  and  compositions  will  continuously  shift 
to  the  equilibrium  condition  at  the  local  pressure.  To 
treat  this  chemistry  problem  rigorously,  one  would  need  to 
incorporate  an  equilibrium  chemistry  computation  into  the 
flow  field  solution.  Rather  than  resort  to  this  complica^ 
tlon,  the  mixing  zone  chemistry  is  treated  by  an  approximate 
but  acceptably  accurate  technique. 

The  equilibrium  mixing  zone  temperatures  and  conw 
positions  are  determined,  as  a  function  of  C,  with  an 
exterior  calculation  for  representative  values  of  p_,  u. 

W  J 

and  u„*.  A  hypothetical  mixing  zone  total  temperature  is 
then  defined,  and  it  is  assumed  that  this  total  temperature, 
as  well  as  the  corresponding  gas  properties,  is  invariant 
throughout  the  flow  field.  This  technique  implies  that  the 
gases  are  in  chemical  equilibrium  at  the  representative 
values  of  p^,  u^  and  u^^,  but  that  the  flow  is  effectively 
frozen  for  deviations  from  the  selected  representative 
condition.  The  validity  of  this  technique  has  been  checked 
for  the  conditions  of  the  rocket-air  mixing  experiments  of 
Refs.  6  and  7,  and  the  error  in  the  mixing  zone  density 
profile  was  found  to  be  negligible. 

Once  the  jet  and  secondary  flow  stagnation  con¬ 
ditions  are  specified,  the  choice  of  a  representative 
value  for  p^  completely  defines  the  static  enthalpy,  H, 

*  The  assumption  of  an  equilibrium  primary  stream  expansion 
process  is  used  to  relate  u,  to  p^. 
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from  eq.  (25a)  and  the  energy  equation: 

A/  =  ”  ^^/z.  (31) 

H,  C  and  completely  define  the  chemical  problem  if  the 
elemental  compositions  of  the  inviscid  streams  are  known. 

The  temperature  and  molecular  composition  can  then  be 
computed  with  a  standard  equilibrium  chemistry  computer 
program.  For  hydrogen— oxygen— nitrogen  or  hydrocarbon— 
oxygen— nitrogen  systems,  the  technique  presented  in  the 
Appendix  of  Ref.  6  may  also  be  used  to  compute  the 
equilibrium  temperatures  and  composition.  A  calculation 
procedure  was  presented  in  Ref.  6  with  which  the  tabulated 
equilibrium  properties  for  hydrogen— air  and  hydrocarbon— air 
systems  (Refs.  37  and  38)  can  be  used  to  compute  the 
equilibrium  properties  for  systems  in  which  the  nitrogen- 
oxygen  mole  ratio  differs  from  3.76  (mole  ratio  for  standard 
air).  That  technique  of  Ref.  6  was  used  in  this  investiga¬ 
tion. 

The  equilibrium  chemistry  computations  yield 
the  static  temperature,  the  molecular  composition  and  the 
molecular  weight  as  a  function  of  The  local  gas  constant, 

R,  is  determined  from  the  molecular  weight,  and  the  mixture 
specific  heat  is  given  by  the  relation 

Cp=^'^I^Cp|^  (32) 

where  m^^  is  the  mass  fraction  of  molecular  species  k  and 
Cp^  is  the  specific  heat  of  species  k  at  the  static  tempera¬ 
ture.  The  specific  heats,  are  tabulated  in  standard 

gas  tables  (Ref.  39,  for  example).  The  local  velocity,  q, 
is  determined  by  U  (from  eq.  25a  and  the  representative 
values  of  u^^  and  Uj),  and  the  hypothetical  mixing  zone  total 
temperature  is  defined  by  use  of  the  energy  equation: 

Ti  =T'  -f-  (33) 
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A  typical  set  of  curves  relating  T^,  c^  and  R  to 

C  is  shown  in  Fig.  4  for  the  experimental  conditions  of 

Refs.  6  and  7.  The  corresponding  frozen  distributions  of 

T  ,  c  and  R  are  also  shown  in  Fig.  4. 
o’  p 

5.  Frozen  Mixing  Zone  Chemistry 

For  frozen  mixing,  two  methods  of  determining  the 

T  ,  c  and  R  distributions  are  presented.  The  first  method 
o  p 

is  used  when  the  specific  heats  of  the  constituent  gases 
may  be  assumed  to  be  constant  and  Independent  of  temperature. 
The  second  method  is  used  when  the  specific  heats  of  the 
gases  are  considered  to  vary  with  temperature.  When  the 
temperatures  of  the  inviscid  streams  are  greatly  different 
the  second  method  should  be  used,  but  the  first  is  adequate 
for  mixing  of  streams  having  relatively  low  temperatures, 
such  as  mixing  in  the  usual  air— air  ejector. 

Constant  specific  heats  —  The  mixing  zone  total 

enthalpy,  H^,  is  determined  from  eq.  (25a)  as  a  function  of 

Cj  where  H _  =  c  T  and  H  .  =  c_.  T  ..  The  local  specific 

oa  pa  oa  oj  pj  oj 

heat  is  given  by 

Cp  —  c")  Cpo, 

The  local  gas  constant  is  given  by 

R  =  CRj 

The  local  mixing  zone  stagnation  temperature  is  given  by 
To  =  Hc/cp 

Variable  specific  heats  —  The  local  mixing  zone 
molecular  composition  is  determined  as  a  function  of  C 
from  the  relation 


(34) 

(35) 

(36) 

(37) 
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Fig.  4  Mixrng  Zone  Properties  for  Rocket-Air  Mixing 
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Where  Is  the  mass  fraction  of  molecular  species  k.  The 
local  gas  constant,  R,  is  given  by  eq.  (35),  but  an  iterative 
process  must  be  used  to  determine  and  Cp.  The  stagnation 
enthalpy,  is  given  by  eq.  (25a),  and  the  static  enthalpy, 

is  computed  from  eq.  (31),  along  with  a  representative  u  vs. 

U  function  (same  procedure  as  for  equilibrium  chemistry) . 

The  local  static  temperature  is  determined  by  trial  of 
various  values  of  T  until  the  following  equation  is  satisfied; 

H  =  ^  (38) 

The  values  of  are  obtained  from  standard  tabulated  gas 
tables  (Ref.  39,  for  example).  After  T  is  determined  the 
mixture  specific  heat  is  obtained  from  eq.  (32) ,  with  the 
c  .  ’s  taken  from  the  gas  tables  at  T.  The  hypothetical 

pK 

total  temperature  is  then  determined  from  eq.  (33) . 

The  resulting  T^,  c„  and  R  functions  of  "C  are 

o  p 

assumed  to  be  invariant  in  the  flow  field  computation.  This 
procedure,  which  is  basically  the  same  as  used  for  equilibrium 
chemistry,  assumes  that  the  variation  of  the  specific  heats 
is  negligible  for  the  range  of  temperatures  that  a  given 
mixture  undergoes  in  the  flow  field.  This  procedure  is 
exact  for  isobaric  mixing,  as  is  true  for  equilibrium 
chemistry. 

The  frozen  distributions  shown  in  Fig.  4  were 
computed  by  taking  into  account  the  variation  of  the 
specific  heats  with  temperature. 
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SECTION  VI 

TURBULENT  EDDY  VISCOSITY 


As  is  usual  with  the  integral  method,  the  turbulent 

shear  is  evaluated  at  only  one  position  in  the  mixing  zone, 

in  this  case  at  the  half  radius  control  surface,  r  (Fig.  3). 

m 

This  turbulent  shear  stress,  t  ,  is  determined  by  the 

m 

equation 

=  fm  1  <®®> 

oV-  f  w\ 

where  is  the  eddy  viscosity. 

Eddy  viscosity  model  —  The  semi— empirical  model  for 
the  eddy  viscosity  which  is  used  in  this  investigation  Is 
the  Prandtl  incompressible  model 

k  b  (40) 

with  the  empirical  constant,  k,  corrected  for  the  influence 

of  variable  density  in  the  mixing  layer.  The  variable 

density  correction  used  Is  a  modification  of  the  correction 

suggested  by  Donaldson  and  Gray  (Ref.  28)  for  the  mixing 

of  an  axlsymmetrlc  jet  flowing  into  a  quiescent  medium. 

Their  analytical  model  was  an  integral  method  similar  to 

the  model  reported  herein,  and  the  turbulent  shear  stress 

was  evaluated  along  a  control  surface,  r^^^  (x) ,  where  the 

velocity  was  u  =  i  u 

^  m  *  max 

Donaldson  and  Gray  evaluated  the  value  of  k  (using 
the  Prandtl  eddy  viscosity  model)  at  various  axial  stations 
in  the  flow  field  for  several  different  experiments.  They 
found  that  the  influence  of  density  gradients  in  the  mixing 
layer  could  be  generalized  if  k  is  taken  to  be  a  function  of 
only  the  local  Mach  number,  M^,  at  the  control  surface,  r^. 

For.  their  cases  of  Interest,  with  zero  secondary  velocity, 
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they  found  an  influence  of  the  nozzle  lip  geometry  on  the 
computed  values  of  k.  For  flows  which  issue  from  a  nozzle 
in  a  large  plate  (causing  the  induced  secondary  flow  to  be 
initially  radial)  the  computed  values  of  k  were,  in  general, 
somewhat  larger  than  for  nozzle  geometries  in  which  the 
induced  secondary  flow  was  nearly  axial.  When  plotted 
against  the  values  of  k  for  both  configurations  tend 
to  decrease  quite  quickly  with  increasing  M  ,  then  to 
become  nearly  constant  for  larger  values  of  M^. 

The  two  curves  which  Donaldson  and  Gray  faired 

through  the  data  (for  sharp  and  blunt  nozzle  lips)  can  be 

approximated  quite  well  by  an  exponential  function  of  M  , 

m 

with  the  ’’half— decay”  of  k  occurring  at  ~  0.45.  If  k^ 
is  the  value  at  •=  0,  and  k  the  asymptotic  value  at 
large  this  exponential  function  can  be  written  as 

=  exp  ^3.42  M>^) 

leo-le« 

Donaldson  and  Gray  evaluated  k  in  both  the  core 
region  of  the  jet  (first  regime)  and  the  developed  region 
(second  regime) .  The  absolute  values  which  they  computed 
for  k  in  the  core  region  appear  to  be  in  error.  Experience 
in  correlating  a  sinillar  integral  method  with  low  speed 
experiments  has  shown  that  the  core  region  (first  regime) 
values  of  k  are  smaller  than  the  developed  region  values  by 
a  factor  of  1.5  —  1.6.  Donaldson  and  Gray,  on  the  other  hand, 
found  the  same  values  of  k  in  both  regions  of  the  flow. 

Their  developed  region  values  of  k  appear  to  be  approxi¬ 
mately  correct,  so  their  core  region  values  are  too  high 


(41) 
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by  a  factor  of  1.5  —  1.6  .  The  reason  for  this  discrepancy 
is  not  clear,  although  It  may  be  related  to  the  momentum 
correction  which  Donaldson  and  Gray  used  to  account  for 
initial  boundary  layers  in  the  jet. 

The  procedure  which  has  been  used  in  this  investi¬ 
gation  is  to  assume  that  eq.  (41)  correctly  expresses  the 
effect  of  density  gradients,  but  to  evaluate  and  by 
correlating  the  present  integral  theory  with  experimental 
results.  The  1-D  Core  Theory  was  correlated  with  various 
low  speed  experiments  (Section  IX)  and  the  values  of  k^  in 
the  first  and  second  regimes  were  found  to  be  approximately 
0.007  and  0.011,  respectively. 

The  2— D  Core  Theory  (Part  II)  was  correlated  with 
several  experiments  for  rocket— air  mixing  (Refs.  6  and  7) 
in  which  varied  from  approximately  1.5  to  1.8.  These 

values  of  are  higher  than  any  reported  by  Donaldson  and 
Gray,  so  the  results  for  k  are  taken  to  be  k  .  It  was  found 
that  the  rocket— air  mixing  could  be  satisfactorily  correlated 
if  kgg/k^  was  taken  to  be  0.66.  Substituting  this  value  into 
eq.  (41) ,  one  obtains 

=  0.66  -f-  O.  34*  exp  (^3,4^ 

This  equation  is  shown  in  Fig.  5,  along  with  the  results 
from  Warren's  experiments,  the  only  experiments  with  sharp 
lip  nozzles  which  Donaldson  and  Gray  considered. 


Because  of  different  definitions  for  the  characteristic 
length  and  velocity  scales  in  the  Prandtl  eddy  viscosity  model 
Donaldson  and  Gray's  values  of  k  are  about  4  times  larger  than 
the  corresponding  values  used  in  this  investigation. 


(42) 
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All  Of  the  computations  presented  in  this  report 
were  made  using  eq.  (42),  along  with  =  0.007  (first 
regime)  and  jjj  =  0.011  (second  and  third  regimes). 

As  was  pointed  out  in  Section  II,  the  Prandtl  eddy 

viscosity  model  (eq.  40)  is  questionable  for  two  stream  mixing 

when  u  /u.  is  greater  than  0.3.  The  greatest  interest  in  this 
a  j 

investigation  has  been  in  correlating  the  rocket-air  mixing 

experiments  of  Refs.  6  and  7.  In  those  experiments,  the 

value  of  u  /u .  did  not  exceed  0.1.  The  2-D  Core  Theory  is 
a  J 

also  correlated  with  supersonic  air-air  ejector  experiments 
in  Part  II,  and  u  /u .  exceeded  0.5  in  some  of  those  experi- 
ments.  It  is  to  be  expected  that  the  predicted  mixing  zone 
growth  rate  will  be  too  small  for  these  air-air  ejectors, 
however,  the  influence  of  mixing  on  the  performance  of 
supersonic  air-air  ejectors  is  known  to  be  small  (in  the 
upstream  choking  mode  of  operation).  Consequently,  the 
use  of  the  Prandtl  eddy  viscosity  model  should  not  cause 
serious  errors  in  the  overall  results. 


Turbulent  shear  stress  —  Differentiating  the 

velocity  profile  equations  from  Section  V,  the  velocity 

gradient  at  r  is  found  to  be: 
m 


First  regime 

_  X  ~  '^S) 

YY)  ^  h 


(43a) 


Second  regime 

ai4|  _  ir  ('u.a.  - 

xr\  t> 


Third  regime 


ax- 


jr, 

a  rw 


(43b) 


(43c) 


Substituting  eqs.  (40)  and  (43)  into  eq.  (39),  one  obtains 
the  following  equations  for  the  shear  stress  at  r^^: 
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First  regime 

Second  regime 

^c<)| 

Third  regime 


U 


^  CW-VJ  “ 


(44a) 


(44b) 


(44c) 


54 


AEOC-TR-68-270 


SECTION  VII 

TRANSFORMATION  OF  DIFFERENTIAL-INTEGRAL  EQUATIONS 


A  system  of  three  basic  differential— integral  equa¬ 
tions  was  derived  in  Section  IV  for  each  of  the  three  mixing 
regimes.  In  this  section,  the  system  of  differential- 
integral  equations  will  be  transformed  into  a  form  suitable 
for  solution  with  conventional  numerical  techniques.  The 
transformation  technique  used  is  essentially  the  same  as  that 
suggested  by  Phares  and  Loper  (Ref.  40). 

1.  First  Regime  Transformation 


The  three  first  regime  differential— integral 
equations  which  were  derived  in  Section  IV  are  as  follows: 


Continuity  Equation 


O  (13) 


Total  Momentum  Equation 

;ri+b 


z  fw -(rj-tt) 


Half— Radius  Momentum  Equation 


■  1=  -  ^ 

J  dX  Z. 


(14) 


Sufficient  information  about  the  terms  in  eqs.  (13) 
through  (15)  was  developed  in  Sections  V  and  VI  so  that  the 


Va') 


(15) 

Z 
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mixing  problem  can  be  defined  in  terms  of  three  flow  field 
variables  which  are  functions  of  only  the  axial  coordinate, 

X.  The  variables  chosen  in  the  first  regime  are  (1)  the  wall 
pressure,  p  ,  (2)  the  radius  of  the  inner  boundary  of  the 

W 

mixing  zone,  r^,  and  (3)  the  width  of  the  mixing  zone,  b. 

The  objective  here  is  to  transform  eqs.  (13)  through  (15) 
into  three  equations  of  the  form 


(45) 


which  may  easily  be  solved  by  standard  numerical  procedures. 

First  apply  Leibnitz*  rule  to  eqs.  (13)  through  (15): 
The  continuity  equation  becomes 


The  total  momentum  equation  becomes 


'T" 

'  c 


The  half-radius  momentum  equation  becomes 


ci?6  Z 

(47) 


+  (vim  -  Uj)  Pj  U.j  'rc‘^^=fr^  (r;+  b^)  _  ^  (r (48) 


Y- 

z 
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The  first  regime  mixing  zone  profiles  are  ex- 

r— rj 

pressed  in  terms  of  the  variable  y  =  — )  (Section  V); 

therefore,  change  the  variable  of  integration  in  the  integral 
terms  from  r  to  y,  using  the  relations 

rl  r 

jf  ')rciT= 


(49) 


r.- 


and 


(  ')rcl'r 


=  (  ' 


d  u  +  b  T, 


(50) 


Equations  (46)  —  (48)  become: 

I 


^  r. 


4- 


j 


-  -  IkL 

dX  ^ 


(47a) 


av. 


=  (r'+b/z)-  C<~t'i>-V2.y 

d  X  2. 


(48a) 
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The  inviscid  stream  parameters  (Pjj  ,  u^) 

are  isentropic  functions  of  only  p  ,  and  the  non-integral 

w 

terms  of  eqs .  (46a)  -  (48a)  can  be  expressed  directly  in 
terms  of  the  flow  field  variables  (x,  p^,  r^,  b) . 

The  integrands  of  the  integral  terms  of  eqs.  (46a)  - 
(48a)  are  S/dx  (pu)  and  d/3x  (pu  ),  and  these  integrands  must 
now  be  expressed  in  terms  of  the  variables  (x,  p^,  r. ,  b) . 

The  parameters  p  and  u  are  actually  functions  of  (x,  r) , 
however,  assume  for  the  moment  that  p  and  u  are  functions 
of  the  independent  variables  (x,  r,  p  ,  r . ,  b) .  Then  one 
can  write 


and 


Where  p  "  "p  and  u  =  u,  the  only  difference  being  the  vari¬ 
ables  which  are  considered  independent.  Application  of  the 
chain  rule  for  differentiation  of  composite  functions  gives 

f  -  ^  4-  ^  ‘=iPvJ  ^  ^  db 


and 


^  - 


dPvj  dX- 

6?^  lijc  dx"  dx 


(51) 


(52) 


Applying  eqs.  (51)  and  (52)  to  the  integrands  6/3x  (pu)  and 
5/3x  (pu^),  one  obtains  _  _ 
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^ (pu,*)  =  lx  tx  ) 


r\  (54) 


Substituting  eqs.  (53)  and  (54)  into  eqs.  (46a) 
through  (48a),  one  obtains  after  expansion  of  the  non-integral 
terms  and  considerable  rearrangement: 


Continuity  equation: 

p  j.  P,  +  r 


(55) 


Total  momentum  equation 

6 


I 


«iPw  ^  ^ 


^  dX 


ci  X 

Half— radius  momentum  equation 


r  r 

SaSx  ^4- 


(56) 


d  X 


(57) 


The  coefficients  of  eqs.  (55)  —  (57)  are  given  by  the 
following  equations:  ^ 


+  j  li-  +  f>.  ^ 

a-  dP./  -“2 


j:  +  u-fc 

^  aiv  i- 


dftt 


['r,N-f<‘c+b)*]^p^ 


(55a) 
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Fa  =  jCf  !?•  +  H- ) d  ^ 

o  o 

F,  =  b^fcr  =^'4-  H 

O  o 

F4 


(55c) 


dry 


w 


w 


d  X 


(55d) 


6,  = 


-  “'■'‘’arw  ^apw 


"IT 


d  iCft, 

d.  pvi/  (56a) 


o  <1  .  (56b) 


o  o  ^ 

G^=  .l,^((2PU.^^^u. 


-  b'C:  ^ 

O 


(56d) 
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H,  =  fpj'i 

°  '4.  °  ^ 

-«^fc"j(f  |^_^+u.  br/j(p^ 


-  +1^^  ^(^M. 

epiw 


^  "^R*^  '■  ■'  ^  a-  515: 


*  y  ^  ' 

-  ««tf l^- +‘^  I?.-)  ^  -  ««  ^’^fc'pf^  .+w|£ ,) 


^  dpw/ 

(57a) 


(57b) 


_  _  '’  ^4. 


*  V  o  ^ 

'Q'jCi® 

o  o'^ 

-htr»(fl-t-  ^^2^  t57d) 
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All  of  the  auxiliary  equations  required  to  evaluate 
the  coefficients  F,  G  and  H  are  presented  In  Appendix  I. 


Equations  (55)  -  (57)  are  linear  in  the  derivatives 
dp^dXj  dr^^/dx  and  db/dx,  therefore  these  equations  can  be 
solved  for  the  derivatives  by  application  of  Cramer's  rule, 
as  long  as  the  determinant  of  the  coefficients  Is  not  zero: 


_  LPj 

dX.  IDI 

-  IR I 
dx  IDI 

d b  ^  \  BI 

ax  lD| 


(58) 

(59) 

(60) 


Where  the  determinants  are: 


Gz  ^3 


F2.  P’3 
Gjy  63-^3 
H4.  Hj.  /V3 


*=^1  P  -4.  F3 

S  \  Ga 


(coefficient  determinant) 


F,  Fz  R4.' 

G,  Gz.  G^\ 

h4 1  \ 


Equations  (58)  —  (60)  can  be  numerically  Integrated  with 
standard  techniques  for  prescribed  initial  conditions. 
The  numerical  technique  used  in  this  investigation  was 
the  well  known  Bunge-Kutta  method  (Ref.  41). 
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2.  Second  Regim*  Transformation 


The  three  second  regime  equations  which  were  de¬ 
rived  in  Section  IV  are  as  follows: 

Continuity  equation 

^  [fu-rc/r  O 

o 


Total  momentum  equation 

o 


Half— radius  momentum  equation 

J.  o 


(16) 


(17) 


(18) 


In  the  second  regime  =  0,  and  the  three  flow 
field  variables  are  chosen  to  be  the  centerline  velocity 
u  ,  and  b.  The  mixing  zone  coordinate  in  the  second  regime 
is  y  “  TT)  therefore 

fc  K-dr  = 

o 


Equations  (16)  —  (18)  are  transformed  in  the  same  way  as 

the  first  regime  equations,  with  the  exceptions  that  the 

variable  u  is  substituted  for  r.,  and  eqs.  (61)  and  (62) 

are  used  in  place  of  eqs.  (49)  and  (50).  One  obtains  the 

following  system  of  equations  for  the  derivatives  dp^dx, 

du  /dx,  db/dx: 
c 


(61) 


(62) 
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Continuity  equation 


<=‘Pw 


ctb  — 

3  air  " 


z.  ^  +p;  =  F. 


■  dX  “ dx 

Total  momentum  equation 

-  dPw  +  5^^  = 


'  dlt  <dx 

Half— radius  momentum  equation 


dK 


w  dPgj  U 

Hi  —I —  -r  rr2-  . 

dX  dX 


The  coefficients  of  eqs.  (63)  —  (65)  are  given  by  the 
following  equations: 


=iPw 


o 

F^.  =-t’‘|CP^+«-3l)v^dv^-R,u.^>r>. 

O 

G,  =  +r^ 


% 


(63) 

(64) 

(65) 

(63a) 

(63b) 

(63c) 

(63d) 

d  Pg, 

dpv^ 

(64  a) 
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|£c')^'i‘^^ 

+u-i£ 

'A.  « 


'6k^C 


Nt-  />' 


■Kzpu. 


db^“-  dbH'^'^r 

- 


'V 


+u^b"[Cp''^+‘^||)>^d^  +r« 


b 

1 


C64b) 


(64c) 


(64d) 


(65a) 


(65b) 


(65c) 


(65d) 
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All  of  the  auxiliary  equations  required  to  evaluate  the 
coefficients  F,  G  and  H  are  presented  in  Appendix  1. 


Equations  (63)  —  (65)  can  be  solved  for  the 
derivatives  of  dp^dx,  du^/dx  and  db/dx  in  the  same  way 
as  in  the  first  regime,  except  that  eq.  (59)  is  replaced  by 


where 


Uil  _  UM 

d  1D| 


Pi  P4  Fs 

Sp 

Hp  H4*  H3 


(66) 


3.  Third  Regime  Transformation 


The  three  third  regime  integral  equations  which 
were  derived  in  Section  IV  are  as  follows; 

Continuity  equation 
■'Tvs/ 


M' 


pdtrcir  =  0 


Total  momentum  equation 


- 

d  X 


XV 

z 


Half— radius  momentum  equation 

O  o  & 


(19) 


(20) 


(21) 


In  the  third  regime  b  =  r^,  and  the  three  flow 

variables  are  chosen  to  be  (p  ,  u„  and  u  ) .  The  mixing 

w  c  w 

zone  coordinate  is  y  =  — ,  therefore 

^w 
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(67) 


and 


r  Vz, 

C  'ird'T  C 

O 


(68) 


Equations  (19)  —  (21)  are  transformed  in  the  same  way  as  the 

first  regime  equations,  with  the  exceptions that  the  variables 

r.  and  b  are  replaced  by  the  centerline  velocity,  u  ,  and  the 
X  c 

wall  velocity,  u^;  eqs.  (67)  and  (68)  are  used  instead  of 
eqs.  (49)  and  (50) .  One  obtains  the  following  system  of 
equations  for  the  derivatives  dp^dx,  du^/dx  and  du^dx; 


Continuity  equation 

^  j,  r*  j-  E"  ^ 


(69) 


Total  momentum  equation 


^  ^  ^  d  LLc  _L  d  iivJ 


die  dDC 

Half— radius  momentum  equation 


4.  4/ 


■f  //j  "il^i  =  i-i. 

d  1C 


The  coefficients  of  eqs.  (69)  —  (71)  are  given  by  the 
following  equations: 

F,  -  ’T.o  ^  f  ‘H- 


r^  =  r\ 


\M 


Tc?  — 


UL 


iL 

^  CCo 


(70) 


(71) 


(69a) 

(69b) 
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_  - 

r-  —  'T  ^  ifr  Na.AH. 

o 

|C^  ^X, 

S,  =  -•-  X 

O 

S^  =  T,m^ 

o 

O 

/4  =  VH 

°  - 


(69c) 


(69d) 


<70a) 


(70b) 


(70c) 


(70d) 


(71a) 


(71b) 
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o 

r''^  -  - 

H^=  -y^^lCs-f^  ^  +‘^’' 

f''^-  - 

+u„,r,j^jCf  ^+  ic  |£  4-tr„  ^w/g, 

o 


All  ol  the  auxiliary  functions  required  to  evaluate  the 
coefficients  F,  G  and  H  are  presented  in  Appendix  I. 


Equations  (69)  —  (71)  can  be  solved  for  the 
derivatives  dp^dx,  du^/dx  and  du^dx  in  the  same  way  as  in 
the  first  regime,  except  that  eq.  (59)  is  replaced  by  eq. 
(66),  and  eq.  (60)  is  replaced  by 


where 


<d  U.W  _  (Jw 


D 


Uw 


F, 

FI  F4. 

*S-| 

H,  Wz.  HJ 


(71c) 


(71d) 


(72) 
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SECTION  VIII 

NUMERICAL  CONSIDERATIONS 


The  system  of  linear  differential  equations  de¬ 
veloped  in  Section  VII,  with  the  coefficients  evaluated  using 
the  equations  of  Appendix  I ,  was  solved  by  conventional 
numerical  procedures.  The  computer  program  is  quite  lengthy 
and  the  details  will  not  be  given  here.  The  purpose  of  this 
section  is  to  give  a  brief  discussion  of  the  numerical  solu¬ 
tion  from  the  point  of  view  of  the  engineering  results. 

Numerical  precision  -  The  integral  terms  appearing 
in  the  coefficients  F,  G  and  H  (Section  VII)  are  evaluated 
by  means  of  a  sixteen  point  Gaussian  quadrature  formula.  No 
improvement  in  the  solutions  was  found  when  the  integrations 
were  made  with  a  32  point  quadrature  formula. 

The  axial  step  size,  Ax,  for  the  Runge-Kutta 
integrations  of  the  set  of  differential  equations  is  an 
input  to  the  program.  The  solutions  have  been  found  insen¬ 
sitive  to  Ax,  and  the  step  size  is  chosen  primarily  for 
convenience  in  defining  the  solution.  Typically,  Ax  is 
chosen  to  be  0.3  —  1.0  times  the  primary  nozzle  radius,  r^. 

Computation  time  —  The  numerical  solutions  were 
obtained  by  use  of  an  IBM— 360/50  digital  computer.  A  typical 
problem  is  solved  in  30  to  40  axial  steps,  requiring  a  cort- 
putation  time  of  approximately  2  minutes. 

Program  inputs  and  printout  —  The  input  parameters 
and  parameters  which  are  printed  out  are  given  in  Appendix 
III. 
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SECTION  IX 

CORRELATION  OF  1-D  CORE  THEORY  WITH  EXPERIMENTS 


•'Most  of  the  correlations  between  theory  and  ex¬ 
periment  have  been  made  with  the  2— D  Core  Theory  (Part  II). 
The  1— D  Core  Theory  has  been  correlated  with  low  speed  air- 
air  experiments  in  order  to  establish  the  incompressible 
values  of  the  eddy  viscosity  constants  (Section  VI)  and  to 
establish  the  validity  of  the  theory  for  constant  density 
flows . 

I.  Low  Speed  Isoboric  Mixing  Experitnentn 

Forstall  and  Shapiro  (Ref>  18)  made  a  series  of 
coaxial  mixing  experiments  with  •=  0.2,  0.25  and  0.5. 

The  experimental  configuration  was  ducted,  but  the  duct 
was  sufficiently  large  in  diameter  so  that  the  induced  axial 
pressure  gradients  were  negligible.  Only  the  first  and 
second  mixing  regimes  were  observed  in  these  experiments. 

The  experimental  lengths  of  the  Inviscid  core  (length  of 
the  first  regime)  are  shown  in  Fig.  6  along  with  the  core 
lengths  predicted  by  the  1-D  Core  Theory  with  =  0.007. 
Also  shown  in  Fig.  6  is  the  experimental  inviscid  core 
length  from  the  zero  secondary  velocity  experiments  of 
Corrsin  and  Uberoi  (Ref.  42).  The  experimental  inviscid 
core  lengths  were  obtained  by  using  the  accepted  procedure 
of  plotting  the  second  regime  centerline  velocity  ratio, 

(u  —  u  y(u.  —  u  )  vs.  X  on  a  log— log  plot  and  then 

C  a  j  a 

extrapolating  the  distribution  back  to  the  axial  position 
where  (u  -  u^)/(u.  —  u  )  equals  one.  Referring  to  Fig.  6, 

C  cl  j  a 

it  is  seen  that  the  theoretical  core  length  is  smaller  than 
the  experimental  value  at  low  secondary  velocity,  but  that 
the  agreement  is  satisfactory  at  velocity  ratios  of  0.2  and 
0.25.  The  theoretical  core  length  is  approximately  25  per 
cent  larger  than  the  experimental  value  at  u^/u.  =  0.5. 

»  j 
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As  discussed  in  Sections  II  and  VI,  this  deviation  from 
experiment  at  velocity  ratios  larger  than  about  0.3  is  normal 
for  the  Prandtl  eddy  viscosity  model. 

The  axial  distribution  of  the  centerline  velocity 

ratio,  (u  -  u^)/(u.  -  u^) ,  is  shown  in  Figs.  7a  —  d  for  the 
c  a  j  a 

Corrsin  and  Uberoi  experiment  and  for  the  experiments  of 

Forstall  and  Shapiro.  The  theoretical  distribution  shown 

in  Fig.  7a  was  computed  at  u^/u.  =  0.01  because  the  theory, 

»  j 

as  currently  formulated,  cannot  be  used  when  the  secondary 

velocity  is  zero*.  Referring  to  Figs.  7a  —  d,  the  theory 

with  =  0.011  predicts  a  somewhat  excessive  rate  of 

decay  for  the  centerline  velocity  ratio  for  u  /u .  %  0,  but 

J 

the  agreement  is  better  for  velocity  ratios  of  0.2  and  0.25. 

As  was  the  case  in  the  first  regime,  the  predicted  rate  of 

mixing  in  the  second  regime  for  u/u.  =  0.5  (Fig.  7d)  is 

a  J 

considerably  less  than  shown  by  the  experiment. 

One  other  aspect  of  the  value  for  k^j  should  be 
noted.  In  the  first  regime  the  mixing  is  initially  two- 
dimensional  in  nature.  One  would  expect,  therefore,  that 
should  be  approximately  the  same  value  as  obtained 
from  two-dimensional  experiments.  Experimental  two- 
dimensional  mixing  data  for  zero  secondary  velocity  are 
usually  correlated  with  the  error  function  profile  (Ref.  26), 

-  '/z  [l  + 


*  The  Corrsin  and  Uberoi  experiment  was  made  without  a  duct. 
The  theoretical  curve  of  Fig.  7a  was  computed  by  assuming 
that  the  process  was  ducted,  but  the  duct  diameter  chosen 
for  the  computation  was  large  enough  so  that  the  resulting 
axial  pressure  gradients  were  negligible . 
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Fig.  7  Continued 
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‘1-D  Core  Theory 

kol  =  0.007,  kQjj  =  0.011 


O  Experiment  (Ref ,  18) 

Forstall  and  Shapiro  Run  C 


( 


40  60 
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where  the  similarity  constant,  a,  is  inversely  proportional 
to  the  rate  of  growth  of  the  mixing  layer.  The  accepted 
incompressible  value  for  cr  is  12.  A  certain  value  for  a 
also  implies  a  certain  value  for  the  rate  of  change  with  x 
of  the  velocity  gradient  at  the  half-velocity  position  in 
the  mixing  zone.  Computations  have  shown  that  the  same  rate 
of  change  of  the  velocity  gradient  at  the  half— velocity 
position  is  obtained  for  the  error  function  profile  with 
a  =  12  and  for  the  cosine  profile  with  k^j  ~  0.007.  This 
result,  therefore,  tends  to  substantiate  the  value  of  k^^ 
which  was  used  to  correlate  the  1— D  Core  Theory  with 
coaxial  mixing  experiments . 

2.  Low  Speed  Jet  Pump  Experiment 

Mikhail  (Ref.  31)  reported  experiments  on  a  low 
speed  jet  pump  with  a  cylindrical  mixing  duct.  The  apparatus 
is  shown  schematically  in  Fig.  8  along  with  the  initial 
conditions  for  the  case  considered  here.  Both  the  primary 
and  secondary  streams  were  unheated  air. 

The  experimental  duct  pressure  distribution  is 
shown  in  Fig.  9  along  with  the  results  of  the  1— D  Core 
Theory.  All  three  mixing  regimes  are  encountered  in  this 
configuration.  The  influence  of  the  wall  boundary  layer  on 
the  axial  pressure  distribution  is  considerable  in  this  long 
mixing  duct.  An  approximate  linear  variation  of  the  wall 
boundary  layer  displacement  thickness,  6  ,  was  used  to 
correct  the  duct  contour  for  the  effect  of  the  boundary 
layer.  This  wall  boundary  layer  correction,  although 
quite  approximate,  causes  the  theoretical  pressure  distri¬ 
bution  to  agree  satisfactorily  with  the  experiment. 

The  centerline  velocity  distribution  for  the  jet 
pump  is  shown  in  Fig.  10.  The  wall  boundary  correction  is 
shown  to  have  a  relatively  small  influence  on  the  center— 
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Fig.  9  Wall  Pressure  Distribution  in  a  Low  Speed  Jet  Pump 
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Centerline  Velocity  Ratio, 


Fig.  10  Centerline  Velocity  Distribution  in  a  Low  Speed  Jet  Pump 
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line  velocity  distribution.  In  these  correlations  the  third 
regime  eddy  viscosity  constant  was  assumed  to  be  the  same 
as  the  second  regime  constant.  The  agreement  between  theory 
and  experiment  indicates  that  the  theoretical  model  for  the 
third  regime  is  reasonably  valid  for  an  axial  distance  into 
that  regime  of  at  least  three  duct  diameters. 

Based  on  these  correlations  of  low  speed  mixing 
experiments,  it  appears  that  the  theory  as  formulated 
provides  a  reasonably  accurate  representation  of  the  con¬ 
stant  density  mixing  process.  It  should  be  noted,  however, 
that  the  range  of  applicability  of  the  Prandtl  eddy  viscosity 
model  is  limited  to  flows  with  low  secondary  velocities 
relative  to  the  primary  velocity. 
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PART  II 

DUCTED  MIXING  THEORY  WITH  TWO-DIMENSIONAL  INVISCID  CORE  FLOW 

(2-D  CORE  THEORY) 

SECTION  X 

DEVELOPMENT  OF  BASIC  DIFFERENTIAL-INTEGRAL  EQUATIONS 

In  Part  I  of  this  report  the  first  regime  flow  was 
solved  by  assuming  that  the  inviscid  core  flow  is  one- 
dimensional.  As  pointed  out  in  Section  I,  the  one— dimensional 
core  flow  assumption  is  seriously  in  error  for  many  practical 
ducted  mixing  problems. 

In  this  sectionj  the  2— D  Core  Theory  differential- 
integral  equations  are  derived  for  the  first  mixing  regime 
(Fig.  11) .  The  mixing  zone  and  secondary  flow  are  treated  In 
essentially  the  same  way  as  in  Part  I,  the  major  difference 
being  in  the  method  of  computing  the  inviscid  core  flow.  The 
basic  differential— integral  equations  are  derived  with  the 
inner  mixing  zone  radius,  r^,  as  the  lower  limit  on  the  inte¬ 
grals.  (It  will  be  recalled  that  the  corresponding  lower  limit 
is  zero  in  the  1— D  Core  Theory  integral  equations.)  A  method 
of  characteristics  solution  of  the  inviscid  core  flow  then 
provides  the  boundary  conditions  (at  r^)  for  the  solution 
of  the  differential— integral  equations. 

1.  Integrotien  of  Boundary  Layer  Equations 

The  fundamental  equations  used  in  this  analysis 
are  the  axial  momentum  equation  and  the  continuity  equation 
(Section  IV): 

Momentum  equation 

PM-v  ^  (1) 

Continuity  equation 
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If  eq.  (I)  is  integrated  between  the  arbitrarily  selected  lower 
and  upper  limits  of  r^  and  T2,  the  following  equation  is 
obtained  (Section  IV); 


= 


Z.  '2. 


dx^ 


(5) 


Note  that  the  pressure,  p^,  is  a  function  of  only  the  axial 
coordinate. 


Integrating  the  continuity  equation  between  the 
limits  and  rg  yields  ^ 

•r. 


(73) 


Substituting  eq.  (73)  into  eq.  (5)  to  eliminate  the  P2V2^2^2 
term: 

4-  (ai-ii-i)?,  \T,  r,  -  a.  A-r 

't/ 


=  n'Tx'-r.'T,  - 


(74) 


Now  let  r 
Equation 


=  r 


1 

lation  (74)  Decomes  ^ 


^  and  Tg 
becomes 


=  r^,  and  note  that  “  0- 


ciPvs/ 


(75) 
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Applying  Leibnitz’  rule  to  eq.  (75),  one  obtains 


dVwj 


dr  a-  3- 

^  I  fU-V +•^‘.■<^1  '('i  -  fwU-wVV, 

'1‘c  ^■Cuj 

f-  1^0  ^  jf u.Tdr-  Uyj  r^-  a:  'Q 


’i  ^  v/ 

.  _  ..  . ,/.. 'i'fv  _  -  ~=tpw  /'^^^v■r.•^^ 

- uv,  _  J 


(76) 


Now  consider  the  continuity  equation  (eq.  73). 
Letting  and  rg  =  r^,  one  obtains  after  applying 

Leibnitz*  rule:  ^ 

Pv,  -  P,- V-,-  =  Jifu.'r(ir-P,-I4;r,-'^ 


^iwW>wJ  '^w 


drv 


vu 


dpC 


but 


thus 


0-^  =.  U.W  ifjii 

ax 

.rvia 


'T* 

Equation  (77)  is  the  desired  form  for  the  continuity 
equation. 


(77) 
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Substituting  eq.  (77)  into  eq.  (76),  one  obtains 
eJ 

ci?C 

which  is  the  desired  form  for  the  "total”  momentum  equation. 


i  ^ 

pu-var  =  fc uc  r-  -  u.i 


a  p,^ 


ax  V  ii  / 


To  derive  the  second  momentum  equation,  again 

consider  eq.  (74).  Letting  r,  =  r.  and  r„  =  one  obtains 
^  X  1  a  m 

r  'Tvii 


applying  Leibnitz*  rule  to  eq.  (79)  one  obtains 

d-X  +  P;  U -  Pn^Uivi  VVi 

+Cay^  -  U.t)  ir^  pLUx-ar--  P;U.;  'T^ 

ax  \  2.  ) 


Uvv\  'Tyy^  ' 


ci  ^ws 


or 


10^ 


jc=‘lL*‘rciv- j  f 

r/ 


(78) 


(79) 


(80) 
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which  is  the  desired  form  of  the  half-radius  momentum  equa¬ 
tion. 


2.  Working  Form  of  the  Differenliol-lnlogrol  Equations 

By  expanding  the  integral  terms  and  noting  that 

r  =  r.  +  b/2,  one  obtains  the  following  set  of  differential— 
m  1 

integral  equations. 


Continuity  equation 
(  I  1  cL 


(81) 


Total 


a# 


momentum 


equation 


Pu-Vir + 


Half— radius  momentum  equation 

^  I  2  J 


puL'ri'f 


(.  -“z — J 

=  (iim-  u.  -Pv 


ax  L  ^  J  <®3) 

One  aspect  of  these  differential— integral  equa¬ 
tions  should  be  noted.  In  the  derivation,  the  boundary 
layer  assumptions  were  assumed  to  apply  in  the  (x,  r)  co¬ 
ordinate  system  (dp/3r  =  0,  etc.).  For  systems  in  which 
the  mixing  layer  is  steeply  inclined  relative  to  the  center¬ 
line,  i.e.,  a  highly  underexpanded  primary  nozzle,  it  would 
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be  more  realistic  to  apply  the  boundary  layer  equations  to 
a  coordinate  system  aligned  with  a  local  characteristic  flow 
direction.  Attempts  were  made  to  derive  the  differential- 
integral  equations  in  such  a  shifting  coordinate  system,  but 
the  resulting  system  of  equations  is  prohibitively  compli¬ 
cated.  The  complexity  is  caused  by  the  necessity  for  matching 
the  mixing  layer  region  with  the  adjacent  regions  of  inviscid 
flow.  It  was  decided,  therefore,  to  use  the  differential- 
integral  equations  as  formulated  here,  even  though  the 
results  will  be  strictly  consistent  with  the  boundary  layer 
assumptions  only  for  small  inclinations  of  the  mixing  layer, 
i.e.,  cos  0.  1.  This  restriction  is  not  serious  for  the 

ducted  mixing  systems  considered  in  this  report;  0^  rarely 
exceeded  15®  —  20®.  Nevertheless,  it  is  important  that 
this  restriction  be  clearly  noted  so  that  the  equations  will 
not  be  applied  indiscriminately  to  mixing  systems  with 
large  inclinations  of  the  mixing  layer. 
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SECTION  XI 

MIXING  ZONE  PROFILES 

The  velocity  and  density  profiles  in  the  mixing 

zone  are  computed  with  the  same  equations  which  were 

developed  in  Section  V  for  the  1— D  Core  Theory.  In  the 

velocity  profile  equation  (eq.  28a),  however,  the  velocity 

u.  is  replaced  by  u.,  defined  as 
3  ^ 

u^  =  cos  0^ 

where  V.  is  the  total  velocity  vector  at  r^  (an  isentropic 
function  of  p^)  and  0^^  is  the  flow  angle  at  r^.  The  method 
of  characteristics  solution  of  the  inviscid  core  flow 
provides  the  inner  boundary  angle,  0^. 

The  velocity  at  the  control  surface,  r^^^,  is 
“m  =  i  t“i  *  '■a> 

Because  one  does  not  integrate  across  the  inviscid 
core  flow  in  the  2— D  Core  Theory,  the  species  conservation 
parameter,  Q,  (eq.  27,  Section  V)  is  evaluated  in  the 
following  manner.  All  of  the  species  which  ^  not  originate 
from  the  central  jet  are  evaluated,  then  this  quantity  is 
subtracted  from  the  total  mass  flow  in  the  system: 


(27a) 
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SECTION  XII 

TRANSFORMATION  OF  DIFFERENTIAL-INTEGRAL  EQUATIONS 


A  system  of  three  differential— integral  equations 
was  derived  in  Section  X  for  the  2— D  Core  Theory.  These 
equations  are  similar  to  the  corresponding  equations  for 
the  1— D  Core  Theory,  and  are  transformed  by  the  same  method 
as  used  in  Section  VII. 


The  three  basic  equations  are  (Section  X) : 


Continuity  equation 

d  r 


(81) 


Total  momentum  equation 


) 


dx  ^  2-  / 


(82) 


Half-radius  momentum  equation 


ax 

n*' 


Pu.*-rd.'r-«.«gj^  [fu-rdr  =■  (Vj +•  *54.) 


Cu.„  -U;)(fiiii'Ci  - 


dpw 

dX 

L  2-  J 

(83) 


The  details  of  the  transformation  are  identical 
to  those  of  Section  VII  and  will  not  be  repeated  here. 
Again,  the  variables  (x,  p^,  r^,  b)  are  considered  indepen- 
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dent.  After  the  transformation,  the  equations  have  the  same 
form  as  in  the  1— D  Core  Theory; 


Continuity  equation 


F,  ^  +  Fi  ^  + 
'  ‘  d?L 


dX 

Total  momentum  equation 


CT  —  CT 


d'T: 


dx  ^  dX. 


6*  =  G. 


Half— radius  momentum  equation 

dp, 


14. 


_:tL  -H  4  * 

dx: 


(55) 


(56) 


(57) 


The  coef f icients .F,  G  and  H,  however,  are  dif¬ 
ferent  from  the  corresponding  1— D  Core  Theory  coefficients 
and  are  given  by  the  following  equations. 


0  O 


4 


Z. 


dfV 


o 


•fw. 


+  i,r.-U||+u, 


(84) 


(85) 
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+■ 


cjy-w 

<=i^ 


(87) 


rt'^7r|  ^  j  ii  ^  A?i 


=  ‘■‘■j’ftfu.a  * 


cl  fV/  a  pw 


(88) 


+  Ar, 


(89) 


I 


l.^f 


(90) 
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°  o^ 

+  n-ii;  ir^r--  -fa,Lt^*'n„  ^ 

///  =  ((zfu-^+u.^ 

r^-  - 


5x)a’J- 


o 


r^- 

o 

- 

‘  _ 


(91) 


(92) 


(93) 


w  I'  r»j 

0  ' 

0 

■  pr:+  b/g.)^-  r^-  ^1 
_  >«  — 

”  r'^_  _  ^*4.  _ 

%  O' 

•  )  '4  +'tn'  j(zf  I 

°  r  '''’■  -  -  °  r'^i-  — 

„bj(p  +  (94) 

o'  ^ 

:  -t="jCzf«.^  +W^|f  -tr;Jc«^^4a^^)dLj. 

+a„,br-  k  |S 

\  ^  ,  X  .  #i«  .ak^  ^ 


o  ^ 

0 

-h(u^^ 


li*“S)=l^ 

UwO  "f"  0"t‘  +■  ^ 
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As  in  Part  I,  eqs.  (55)  -  (57)  are  linear  in  the 
derivatives  dp^dx,  dr^/dx  and  db/dx,  and  can  be  solved 
for  the  derivatives  by  use  of  Cramer's  rule  (eqs.  58—60, 
Section  VII).  Equations  (58)  —  (60)  are  again  numerically 
integrated  by  use  of  the  Runge— Kutta  method. 

Many  auxiliary  functions  are  necessary  to  evaluate 
the  coefficients  F,  G,  and  H.  Those  functions  which  are  not 
the  same  as  those  of  the  1— D  Core  Theory  (Appendix  I)  are 
presented  in  Appendix  II.  It  is  important  to  note  that 
the  coefficients  contain  terms  which  depend  on  the  flow 
angle,  9^^,  at  r^.  The  local  value  of  9^^  is  computed  in 
the  method  of  characteristics  solution  of  the  inviscid 
core  flow  (Section  XIII). 
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SECTION  XIII 

METHOD  OF  CHARACTERISTICS  SOLUTION  OF  INVISCID  CORE  FLOW 


Many  terras  appearing  in  the  coefficients  F,  G  and 
H  (eqs.  84—95,  Section  XII)  depend  on  the  flow  conditions 
at  the  inner  mixing  zone  boundary.  In  order  to  solve  the 
system  of  equations,  one  must  be  able  to  evaluate  the 
following  parameters  at  r^:  0^,  du/3x,  du/Sr,  3p/dr,{dp/Bx  (see 

Appendix  II).  To  provide  these  parameters,  the  inviscid 
core  flow  is  developed  (with  the  irrotational  method  of 
characteristics)  simultaneously  with  the  numerical  solution 
of  the  equations  for  the  mixing  layer  and  inviscid  secondary 
flow. 


1.  Interior  Flow  Salufian 

The  method  of  characteristics  is  a  well  established 
method  for  computing  supersonic  inviscid  flow  fields.  De¬ 
tailed  discussions  of  the  method  of  characteristics  are 
presented  in  standard  references  on  supersonic  aerodynamics 
(Ref.  43,  for  example).  Briefly,  the  flow  field  solutions 
are  obtained  by  computing  along  ’’characteristic”  lines  in 
the  flow  field  which  correspond  to  Mach  lines.  The  change  of 
flow  conditions  along  these  characteristic  curves  is  computed 
by  solving  the  following  set  of  ordinary  differential  equa^ 
tions  (for  axisymmetric  Irrotational  flow) : 


—  M  .  Sim  ^ 

<X  4-  ■  ,  _ 


tcur\  (b  +  0^^ 


(97) 
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The  upper  sign  refers  to  right  running  (Family  I)  charac¬ 
teristics  and  the  lower  sign  to  left  running  (Family  II) 
characteristics.  The  local  Mach  angle  is  denoted  by  a. 

Consider  a  unit  process  in  the  physical  plane 
(Fig.  12).  If  the  flow  conditions  are  known  at  points  1  and 
2  (which  are  not  on  the  same  characteristic),  then  eqs.  (96) 
and  (97)  can  be  solved  to  locate  point  3  (Xg,  r^)  and  to 
define  Vg  and  Qg.  For  a  small  unit  process  the  Mach  lines 
are  assumed  to  be  straight  and  eqs.  (96)  and  (97)  are 
written  in  finite  difference  form.  As  a  first  approximation 
one  evaluates  the  coefficients  in  the  difference  equations 
at  the  initial  conditions.  The  unit  process  is  then  re¬ 
calculated  using  average  values  for  the  coefficients  until 
the  solution  converges. 

In  the  inviscid  core  solution,  the  initial  condi¬ 
tions  are  prescribed  along  a  radial  line  AB  at  the  primary 
nozzle  exit  plane  (Fig.  13).  For- an  underexpanded  primary 
nozzle,  the  flow  will  expand  at  the  nozzle  lip  through  a 
Prandtl— Meyer  expansion  until  the  local  pressure  is  equal 
to  If  AC  represents  the  Family  I  characteristic  which 

originates  as  the  last  Mach  line  in  the  Prandtl— Meyer  ex¬ 
pansion  fan,  then  the  flow  in  the  region  ABC  can  be  corti- 
puted,  based  on  only  the  initial  conditions.  Note,  however, 
that  the  extent  of  region  ABC  depends  on  the  value  of  . 

2.  Boundary  Solution 

Consider  the  flow  near  the  inner  mixing  zone 
boundary  at  some  downstream  location,  x  (Fig.  14).  At  x, 
the  upstream  solution  has  defined  the  values  of  p^  and  r^^. 
The  derivatives  dp^dx  and  dr^/dx  are  also  specified  at  x. 
Since  the  flow  is  completely  defined  at  x,  a  new  Family  I 
characteristic  can  be  generated  from  the  boundary  point  D. 
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Fig.  12  Unit  Process  in  Method  of  Characteristics  Solution 
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Family  I  Characteristic 
Originating  as  the  Last 
Line  in  Prandtl-Meyer 
Expansion  Fan 


Fig.  13  Method  of  Characteristics  Solution  of  Inviscid  Core  Flow 
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Fig.  14  Solution  at  o  General  Boundary  Point 
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(The  "interior”  solution  is  thus  extended  outward  by  one 
Family  I  characteristic  at  each  computation  interval,  Ax). 

With  the  conditions  at  x  completely  known,  the 
Runge— Kutta  method  prescribes  a  series  of  intermediate  points 
at  which  dp^dx,  dr^/dx  and  db/dx  are  to  be  computed. 

These  specified  points  (x,  p  ,  r . ,  b)  do  not  necessarily 
correspond  to  points  on  the  actual  solution,  but  the  solu¬ 
tion  depends  on  a  weighted  combination  of  the  intermediate 
results.  Such  a  Runge— Kutta  point  is  shown  as  point  E  in 
Fig.  14.  With  the  position  of  E  specified,  along  with  the 
pressure,  the  only  unknown  at  E  is  the  angle,  6^.  The  angle, 
0^,  which  is  required  for  the  succeeding  Runge— Kutta  points, 
is  determined  by  the  following  procedure.  A  Family  II 
characteristic  EF  is  constructed fWhich  intersects  the 
last  Family  I  characteristic,  along  which  the  flow  is 
completely  defined  from  the  interior  solution.  By  an 
Iterative  solution  of  eqs.  (96)  and  (97),  the  location 
of  F  and  the  value  of  0^  at  E  are  determined. 

After  the  system  of  equations  is  solved  lor  p^ 
and  r^  at  X  +  Ax  (point  G,  Fig-  14),  the  above  procedure 
is  repeated  to  define  the  physically  correct  value  of  9^, 
at  X  +  Ax. 

3.  Inner  Boundary  Derivatives 

The  value  of  0.  at  each  Runge— Kutta  point  is 
computed  with  the  procedure  described  in  the  preceding 
section.  In  addition  to  0^,  the  derivatives  du/3x, 
du/0r,  3p/3x  and  3p/9r  at  r^  must  be  computed  at  each 
Runge— Kutta  point. 

Now  consider  the  following  system  of  equations 
which  defines  the  flow  in  the  inviscid  core  region; 
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JoL  =.  dx 

+  ¥7dr 

(98) 

-1-  ^  <ir* 
dr- 

(99) 

dp- 

4-  ^  dr* 
dr 

(100) 

pui  1^  -1-  Pw- 

dUj  ^ 

^r  - 

(101) 

_ 

dr 

(102) 

|^(f  u.r)  + 

|^(e<rr)=o 

(103) 

After  expansion  and  introduction  of  the  equation  of  state, 
the  continuity  equation  (eq-  103)  may  be  written: 


P0“  (103a) 


The  above  system  of  equations  :±s  linear  in  the  derivatives 

and  is  solved  for  these  derivatives  at  each  Bunge-Kutta 
point  by  means  of  standard  matrix  techniques  . 


This  method  of  computing  the  inner  mixing  zone  boundary 
derivatives  was  suggested  by  the  work  of  Loper  and  Lightsey 
(Ref.  44). 


102 


AEDC-TR-dB.270 


Referring  to  Fig.  14,  the  total  differentials 
are  represented  by  the  following  differences: 

dx  =  xlg-  Xt, 

dr  'fie-'Qo 

d  U.  =  -  til  D 

dss  ^  D 

d  p  =  RwC  -  D 

Treating  these  total  differentials  as  differences  is 
admittedly  approximate;  however,  the  accuracy  of  this 
approximation  is  consistent  with  the  accuracy  of  the 
method  of  characteristics  solution  which  involves  similar 
approximations . 
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SECTION  XIV 

NUMERICAL  CONSIDERATIONS-2.D  CORE  THEORY 

The  computer  program  for  the  2— D  Core  Theory  is 
much  more  complex  than  the  program  for  the  1— D  Core  Theory, 
and  the  details  of  the  program  will  not  be  presented  here. 

The  purpose  of  this  section  is  to  give  a  brief  discussion 
of  the  numerical  solution  from  the  point  of  view  of  the 
engineering  results. 

Axial  Step  Size  —  The  axial  step  size,  Ax,  for  the 
Runge—Kutta  integration  of  the  system  of  differential  equa^ 
tions  must  be  much  smaller  than  the  corresponding  step  size 
in  the  1— D  Core  Theory.  The  step  size  is  chosen  to  be  about 
the  same  size  as  a  typical  grid  in  the  method  of  charac¬ 
teristics  network  (typically  0.05—0.10  times  r^) .  If  the 
step  size  were  to  be  larger,  the  accuracy  of  the  "boundary 
solution"  (Section  XIII)  would  be  poor. 

Method  of  Characteristics  —  The  accuracy  of  the 
method  of  characteristics  is  well  established  as  long  as 
the  network  is  made  reasonably  small  compared  to  the  local 
rate  of  change  of  the  flow  conditions.  One  point  of  interest 
is  that  the  network  of  points  is  redistributed  along  a 
radial  line  (by  interpolation)  at  each  axial  station  in  the 
solution.  By  this  redistribution  procedure,  the  problem  of 
network  "foldback",  indicative  of  the  onset  of  boundary 
shock  formation,  is  avoided.  The  method  of  characteristics 
portion  of  the  program  has  been  compared  with  standard  "jet 
plume"  programs  and  the  results  are  nearly  Identical. 

Computation  Time  —  The  numerical  solutions  were 
obtained  by  use  of  an  IBM— 360/50  digital  computer.  A 
typical  problem  is  solved  in  100  to  120  axial  steps, 
requiring  a  total  time  of  approximately  10  minutes. 
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Program  Inputs  arid  Printout  -  The  input  parameters 
and  parameters  which  are  printed  out  are  given  in  Appendix 
III. 
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SECTION  XV 

REGIMES  OF  SOLUTIONS 


The  method  of  computation  depends  on  the  mode  of 
operation  of  the  ducted  mixing  system  (Fig.  2,  Section  I). 

In  all  cases,  the  procedure  is  to  assume  a  value  of  the 
initial  static  pressure,  »  make  a  flow  field  computation, 
and  then  iterate  on  p^j^  until  some  downstream  criterion  is 
satisfied.  In  the  "downstream  choking"  and  "back  pressure 
dependent"  modes  the  criterion  is  either  that  the  flow  must 
choke  at  the  duct  exit  or  that  the  duct  exit  pressure  must 
match  a  prescribed  back  pressure. 

The  problem  is  much  more  complex  in  the  case  of 
the  "upstream  choking"  mode  (Fig.  2a,  Section  I).  In  this 
mode  the  flow  must  accelerate  to  the  critical  condition,  and 
then  the  solution  must  pass  smoothly  into  the  supercritical 
regime.  When  the  flow  becomes  critical,  the  determinant  |d| 
for  the  system  of  equa'tions  (Section  VII)  is  zero.  The 
basic  premise  which  is  used  to  compute  the  correct  solution 
is  that  the  static  pressure  is  continuous  through  the  criti¬ 
cal  region;  therefore,  one  Iterates  on  p^j^  until  the  deter¬ 
minants  |d|  and  |p[  (Section  VII)  simultaneously  become  zero. 
The  technique  for  finding  the  correct  solution  is  illustrated 
in  Fig.  15.  If  one  chooses  a  value  of  p^j^  which  is  higher 
than  the  correct  value,  the  flow  will  accelerate  and  then 
decelerate,  with  |d|  remaining  positive.  If  the  selected 
p^l  is  too  low,  then  |D[  will  go  to  zero  upstream  of  the 
position  where  |p|  goes  to  zero. 

Although  not  required  by  this  solution  technique, 
the  determinants  |r|  and  |b|  (Section  VII)  also  go  to  zero 
at  the  critical  section.  Thus  the  parameters  r^  and  b  are 
continuous  through  the  critical  region. 
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Axial  Distance,  x 


Fig.  T5  Solutions  for  Ductod  Mixing  System  Oporoting  in  Upstream  Choking  Mode 
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In  principle,  one  needs  only  to  define  the  correct 

value  of  p^j^  within  sufficiently  small  limits  and  then  the 

solution  will  run  smoothly  through  the  critical  point.  By 

luck,  the  first  solution  attempted  ran  smoothly  without  a 

very  small  tolerance  on  It  was  later  found  that,  in 

general,  the  solutions  would  not  run  smoothly  even  with  a 

tolerance  on  p^^  as  small  as  one  part  in  a  million.  To 

eliminate  the  necessity  for  the  extremely  large  number  of 

iterations  required  to  bracket  p^^^  to  the, required  tolerance, 

the  following  procedure  was  used.  After  bracketing  p  ,  within 

4 

a  reasonably  small  range  (to  one  part  in  10  ,  typically) ,  the 
solution  is  rerun  for  at  the  average  value  of  the  last 

’'low"  and  ’’high”  solutions  (Fig.  15).  The  secondary  stream 
total  pressure  is  arbitrarily  increased  by  a  small  amount 
(typically  1%)  at  the  axial  computation  station  immediately 
upstream  of  the  position  where  the  flow  tends  to  become 
critical.  The  solution  then  runs  smoothly  into  the  super¬ 
critical  regime.  This  method  of  computing  the  supercritical 
solutions  has  been  automated  in  the  computer  program. 

If  one  does  not  know  beforehand  in  which  mode  that 
a  particular  ducted  system  will  operate,  the  procedure  is  to 
select  a  high  initial  value  of  Successive  computations 

are  then  made  at  lower  and  lower  values  of  ‘“'^til  (1)  the 

duct  exit  pressure  matches  the  prescribed  back  pressure, 

(2)  the  flow  chokes  at  the  duct  exit  (indicated  by  |d|  =  0) 
or  (3)  the  flow  becomes  critical  at  some  upstream  position 
in  the  duct. 

It  should  be  noted  that  the  above  procedures  are 
also  applicable  to  the  1-D  Core  Theory,  although  the  upstream 
choking  mode  would  not  ordinarily  be  encountered  if  the 
inviscid  core  flow  is  assumed  to  be  one-dimensional. 
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SECTION  XVI 

CORRELATION  OF  2-D  CORE  THEORY  WITH  AIR-AIR  EJECTOR  EXPERIMENTS 


The  2— D  Core  Theory  has  been  correlated  with  the 
air— air  ejector  experiments  of  Chow  and  Addy  (Ref.  5)  and 
of  Chow  and  Yeh  (Ref.  45).  All  of  the  experimental  cases 
considered  here  were  operated  in  the  upstream  choking  mode, 
l.e.,  with  relatively  low  back  pressure.  The  experimental 
configurations  are  shown  in  Fig.  16. 

1.  Constant  Area  Air-Air  Ejectors 

Chow  and  Addy  (Ref.  5)  presented  a  series  of 
experiments  for  the  cylindrical  air— air  ejector  shown  in 
Fig,  16a.  The  primary  nozzle  was  designed  to  provide  a 
uniform  primary  flow  at  =  2.  Only  the  nozzle  configura¬ 
tion  with  the  thinnest  lip  will  be  considered  here.  Both 
the  primary  and  secondary  streams  were  unheated  air. 

The  experimental  mass  flow  ratio,  w^/Wj ,  is  shown 
in  Fig.  17,  along  with  the  results  of  the  2— D  Core  Theory. 

The  agreement  between  theory  and  experiment  is  very  satis¬ 
factory  for  the  entire  range  of  p  over  which  the 

theory  is  applicable.  At  total  pressure  ratios  below  about 
0.05,  the  mixing  zone  impinges  on  the  wall,  and  the  2— D  Core 
Theory  as  presently  formulated  is  not  valid.  To  compute 
these  low  secondary  mass  flows,  Chow  and  Addy  used  the  base 
pressure  analysis  of  Korst,  et  al.  (Ref.  26),  suitably 
modified  to  include  base  bleed  effects.  Such  a  base  pres¬ 
sure  analysis  was  not  incorporated  into  the  2— D  Core  Theory 
because  the  main  interest  is  in  ducted  mixing  with  relatively 
large  secondary  flows.  The  theoretical  results  of  Chow  and 
Addy  are  not  shown  in  Fig.  17;  they  are  nearly  identical  to 
the  2— D  Core  Theory  results. 
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c.  Disploced  Pgrabolic  Duci  of  Chow  and  Yeh  (Ref.  45) 
Fig.  16  Air-Air  Ejector  Configurotions 
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Fig.  17  Mojs  Flew  Ratio  for  a  Cylindrical  Air— Air  Ejector 
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The  initial  secondary  Mach  number,  is  plotted 

against  the  pressure  ratio  P,„t/p^.  in  Fig.  18.  Again,  the 
agreement  between  the  2-D  Core  Theory  and  experiment  is  very 
good  over  the  range  of  applicability  of  the  theory.  For 
p^j^  >  PjjL’  ^al  ~  ^  which  corresponds  to  Fabri  and  Paulon’s 
"saturated  supersonic  regime." 

The  wall  pressure  distributions  from  Chow  and  Addy's 
experiments  are  shown  in  Fig.  19.  The  2-D  Core  Theory  predicts 
the  experimental  distributions  fairly  well  for  all  three  pres¬ 
sure  ratios  shown.  Also  shown  in  Fig.  19  is  the  locus  of 
positions  where  the  theoretical  flow  becomes  critical.  It  is 
noteworthy  that  the  secondary  Mach  number  is  subsonic  at  the 
critical  point ,  and  becomes  more  subsonic  as  the  secondary  flow 
is  decreased.  The  system  of  equations  is  written  for  the 
combined  supersonic  and  subsonic  portions  of  the  flow,  and 
the  results  reflect  that  if  part  of  the  flow  is  supersonic 
at  the  critical  point,  then  another  part  of  the  flow  must  be 
subsonic.  This  phenomenon  is  consistent  with  some  of  the 
recent  results  on  the  choking  of  a  nonuniform  flow  (Ref.  46, 
for  example)*. 

2.  Variable  Area  Air-Air  Ejectors 

Chow  and  Yeh  (Ref .  45)  presented  a  series  of  experi¬ 
ments  for  an  air-air  ejector  in  which  the  wall  shape  was 
parabolic  (Fig.  16b) .  They  also  made  experiments  with  the 
parabolic  duct  displaced  downstream,  so  that  the  minimum  duct 
area  was  located  2r^  downstream  from  the  primary  nozzle  exit 
(Fig.  16c). 

The  mass  flow  ratio,  w^/w.,  is  shown  in  Fig.  20 

^  J 

for  the  parabolic  duct  (Fig.  16b).  As  was  the  case  for  the 
cylindrical  duct,  the  2-D  Core  Theory  predicts  the  experi¬ 
mental  results  quite  accurately  over  the  entire  range  of 
applicability. 

*  The  interpretation  that  )Dj  =  0  corresponds  to  physical  choking 
is  strengthened  by  the  fact  that  |d|  =  0  occurs  almost  exactly  at 
M  =  1  for  coaxial  flows  computed  with  the  assumption  of  zero 
mixing  zone  shear  stress  (no  mixing  case) . 
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Fig.  20  Mast  Flow  Ratio  for  o  Voriablo  Area  Air-Air  Ejector 
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The  initial  secondary  Mach  number,  is  plotted 

against  the  pressure  ratio  p,„n/p_„  in  Fig.  21  for  both  the 
parabolic  duct  and  the  displaced  parabolic  duct.  Because  of 
the  downstream  displacement  of  the  minimum  area,  the  saturated 
supersonic  regime  is  not  encountered  with  the  displaced 
parabolic  duct . 

The  wall  pressure  distributions  for  the  displaced 
parabolic  duct  are  shown  in  Fig.  22.  Some  of  the  experiments 
were  made  with  fairly  high  back  pressure,  and  these  experi¬ 
mental  distributions  show  that  the  back  pressure  has  caused 
the  wall  boundary  layer  to  separate  downstream  of  the  criti¬ 
cal  point  (especially  for  “  0.06).  The  2-D  Core 

Theory  predicts  the  experimental  distributions  fairly  well, 
especially  at  the  lower  values  of  P^2.^Poa  values  of 

V*j>- 

Based  on  these  correlations  of  air-air  ejector-’ 

experiments,  It  appears  that  the  2— D  Core  Theory  accurately 

represents  the  flow  processes  in  such  configurations.  Even 

though  u  /u.  exceeded  0.5  in  the  downstream  portions  of  the 
a  j 

flow  fields,  the  deficiency  in  the  Prandtl  eddy  viscosity 
model  is  not  apparent  from  these  correlations.  As  pointed 
out  by  Chow  and  Addy,  the  flow  in  such  air— air  ejectors 
(with  relatively  large  secondary  mass  flows)  is  predominately 
controlled  by  the  inviscid  stream  Interactions.  The  effects 
of  mixing  are  small  enough  so  that  errors  in  the  mixing  rate 
are  not  readily  discernible. 

Because  the  viscous  effects  are  small  in  these 
supersonic  air— air  ejector  configurations,  the  accurate 
prediction  of  the  experimental  performance  depends  primarily 
on  accurately  computing  the  inviscid  flow  fields.  These 
correlations  of  the  2— D  Core  Theory  with  experiment  indicate, 
therefore,  that  the  present  method  of  computing  the  inviscid 
flowsis  satisfactory,  at  least  in  the  near  field  (x<5— 6  r^) . 
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Fig.  22  Wall  Pressure  Distributions  in  a  Variable  Area  Air— Air  Ejector 
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SECTION  XVII 

CORRELATION  OF  2-D  CORE  THEORY  WITH  ROCKET-AIR  MIXING  EXPERIMENTS 


The  main  objective  of  this  investigation  has  been 
to  predict  the  performance  of  ducted  mixing  systems  in  which 
chemical  reactions  occur.  Peters,  et  al.  (Ref.  6)  and 
Cunningham  and  Peters  (Ref.  7)  presented  experimental  re¬ 
sults  on  ducted  mixing  systems  in  which  the  primary  stream 
was  the  exhaust  from  a  hydrogen^oxygen  rocket  and  the 
secondary  stream  was  unheated  air.  The  rocket  was  operated 
with  approximately  2.5  times  stoichiometric  fuel  (0/F  =  3.2), 
therefore,  there  was  potential  for  considerable  heat  release 
in  the  mixing  process.  The  rocket  nozzle  area  ratio  was 
approximately  5,  and  the  rocket  was  operated  with  a  com¬ 
bustion  chamber  pressure  (pQ.j)  of  approximately  20  atm,  for 
all  tests.  The  two  mixing  duct  configurations  shown  in 
Fig.  23  will  be  considered  here.  Experimental  results  were 
presented  in  Ref.  6  for  the  conical  mixing  duct  of  Fig.  23a. 
The  mixing  system  was  operated  in  the  upstream  choking  mode 
(Fig.  2,  Section  I).  Results  were  presented  in  Ref.  7  for 
the  same  basic  apparatus,  except  that  a  cylindrical  exten¬ 
sion  was  added  to  the  conical  mixing  duct  (Fig.  23b)  so  that 
the  system  operated  in  the  downstream  choking  mode.  Results 
were  also  presented  in  Ref.  7  for  the  configuration  of 
Fig.  23b  operated  in  the  back  pressure  dependent  mode. 

References  6  and  7  present  experimental  wall 
pressure  distributions,  secondary— primary  mass  flow  ratios, 
mixing  duct  thrust  ratios,  and  radial  distributions  of  Pitot 
pressure  and  gas  composition  at  the  duct  exit  plane. 

The  mixing  duct  thrust,  F^,  is  defined  as  the 
axial  pressure  force  exerted  on  the  mixing  duct  between  the 
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r^/r^l  “  0.38,  r^i  =  3  In.  {«  7.5  cm.) 
Conical  Rocket  Nozzle  (15^^  Half  Angle) 


a.  Conical  Mixing  Duct  (Ref.  6) 


b.  Conical  Mixing  Duct  with  Cylindrical  Extension  (Ref.  7) 
Fig.  23  Rocbet-Air  Mixing  Configurations 


AEDC-TR-6a-270 


entrance  and  exit  sections,  i.e., 

r  ^ 

R  =2Tr  (104) 

Where  L  is  the  length  of  the  mixing  duct  and  is  a  known 
function  of  x. 

1.  Mixing  Syslem  Operating  in  th«  Upstrnam  Choking  Mod* 

This  section  presents  the  correlation  of  the  2— D 
Core  Theory  with  the  experimental  results  for  the  configura¬ 
tion  of  Fig.  23a  (Ref.  6). 

Axial  distribution  of  theoretical  flow  parameters  —  Before 
discussing  the  correlations  of  the  2— D  Core  Theory  with 
experiments,  consider  the  theoretical  distributions  of  the 
flow  field  parameters  shown  in  Fig.  24.  The  axial  distribu¬ 
tions  of  p  ,  r . ,  b  and  Q  are  shown  in  Fig.  24a  for  a  typical 

W  ^ 

ducted  mixing  computation.  The  distributions  of  6^  and  the 
coefficient  determinant,  |D|  ,  are  shown  in  Fig.  24b.  This 
flow  field  solution  was  obtained  by  using  the  procedure 
described  in  Section  XV,  i.e.,  the  initial  static  pressure, 

was  iterated  until  the  solution  was  bracketed  to  a  small 
tolerance.  The  secondary  stream  total  pressure,  p__,  was 

Oo. 

then  increased  slightly  (by  0.3%  in  this  case)  to  allow  the 
solution  to  proceed  smoothly  into  the  supercritical  regime. 

The  pressure  distribution,  shown  in 

Fig.  24a  indicates  that  the  secondary  Mach  number,  M„ ,  is 

O' 

initially  0.71;  increases  to  approximately  0.94  at  the 
critical  section.  The  value  of  at  the  duct  exit  is 
approximately  1.3. 

The  axial  distribution  of  (Fig.  24a)  shows  the 
usual  tendency  of  a  jet  from  a  conical  nozzle  to  have  a 
periodic  structure.  If  the  flow  were  inviscid,  the  second 
peak  in  r^  would  be  about  as  large  as  the  first  peak.  With 
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mixing,  the  core  flow  tends  to  expand  and  contract,  but  with 
gradually  decreasing  r^.  This  gradual  decrease  In  r^  is 
caused  by  the  mixing  zone  growing  into  the  inviscid  core  flow. 
If  the  duct  were  sufficiently  long  the  inviscid  flow  would 
eventually  be  completely  dissipated.  With  this  configuration, 
however,  only  the  first  mixing  regime  is  encountered. 

The  axial  distribution  of  b  (Fig.  24a)  shows  that 
the  mixing  zone  growth  rate,  db/dx,  is  not  approximately 
constant,  but  is  significantly  influenced  by  the  periodic 
nature  of  the  inviscid  core  flow.  Computations  have  been 
made  for  similar  pressure  distributions  with  the  1—D  Core 
Theory,  and  no  large  variation  of  db/dx  was  predicted.  It  is 
concluded,  therefore,  that  the  variations  of  db/dx  indicated 
by  Fig.  24a  are  primarily  caused  by  the  angularity  of  the 
flow  at  r^.  There  is  a  definite  qualitative  relation  between 
the  distribution  of  0^  (Fig.  24b)  and  b,  with  negative  values 
of  tending  to  increase  the  mixing  zone  growth  rate. 

The  axial  distribution  of  the  coefficient  deter¬ 
minant,  |D|  ,  for  the  system  of  equations  (Fig.  24b)  shows 
that  |Dj  is  positive  for  subcritlcal  flow  and  passes  smoothly 
through  zero  at  the  critical  point. 

The  species  conservation  parameter,  Q,  (defined  by 
eq.  2V. ,  Section  XI)  is  shown  in  Fig.  24a.  It  will  be  recalled 
that  if  Q  does  not  deviate  from  unity,  then  the  use  of  eq. 

(25a)  (Section  V)  to  relate  the  velocity,  composition  and 
enthalpy  profiles  is  justified.  It  is  seen  that  Q  is 
initially  unity,  then  undergoes  a  step  decrease  to  about  0,98 
at  the  critical  section  where  p^^  was  arbitrarily  increased. 
Downstream  of  the  critical  section,  Q  gradually  increases  and 
has  a  value  of  1.01  at  the  duct  exit-  This  behavior  of  Q  is 
typical  for  all  of  the  rochet-air  mixing  cases  considered. 

In  no  case  does  Q  deviate  from  unity  by  more  than  about  2%. 
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It  can  be  concluded,  therefore,  that  eq.  (25a)  satisfactorily 
relates  the  species  and  enthalpy  profiles  to  the  mixing  zone 
velocity  profiles,  at  least  for  use  in  an  integral  method. 

Wall  pressure  distributions  -  Experimental  wall  pressure 
distributions  are  shown  in  Figs.  25a  and  25b  for  two  values 
of  the  total  pressure  ratio,  Pq^/Pqj •  Also  shown  are  the 
pressure  distributions  predicted  by  the  2— D  Core  Theory  with 
the  assumption  of  equilibrium  mixing  zone  chemistry.  The 
theoretical  pressure  distribution  agrees  fairly  well  with  the 
experiment  in  the  upstream  half  of  the  duct,  but  falls  some¬ 
what  below  the  experiment  in  the  downstream  part  of  the  flow. 
One  should  note  the  tendency  of  the  Inviscid  core  flow  to 
alternately  expand  and  contract  (Fig.  24a)  which  has  a 
corresponding  effect  on  the  wall  pressure  distribution.  The 
theoretical  second  peak  in  r^^  is  apparently  too  small.  This 
deviation  from  experiment  may  be  explained,  at  least  quali¬ 
tatively,  by  the  presence  of  shock  waves  in  the  experimental 
core  flow.  If  the  total  pressure  of  the  core  flow  were  re¬ 
duced,  on  the  average,  by  a  system  of  shock  waves,  then  the 
core  flow  would  occupy  a  larger  area  in  order  to  satisfy 
continuity  relations.  The  greater  area  occupied  by  the  core 
flow  would  decrease  the  area  available  for  the  supersonic 
secondary  flow,  therefore,  the  theoretical  wall  static  pres¬ 
sure  would  be  increased  in  the  downstream  part  of  the  flow. 

A  theoretical  pressure  distribution  is  also  shown 
in  Fig.  25a  for  the  case  of  no  mixing.  Upstream  of  the  criti¬ 
cal  section,  the  wall  pressure  is  slightly  lower  than  in  the 
case  with  mixing.  Because  of  the  higher  the  secondary 

mass  flow  is  approximately  3%  higher  for  the  zero  mixing 
case.  The  theoretical  pressure  distribution  for  zero  mixing 
is  considerably  lower  than  the  theoretical  distribution  with 
mixing  in  the  downstream  part  of  the  flow.  Consequently,  the 
mixing  duct  thrust,  F^,  is  about  25%  lower  than  the  value 
obtained  with  mixing. 
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Comparing  the  two  theoretical  pressure  distributions 
in  Fig.  25a,  one  sees  that  the  mixing  acts  as  blockage  to  the 
flow.  An  effective  displacement  thickness,  6*,  for  the 
mixing  layer  is  then  positive.  It  is  noteworthy  that  for 
mixing  in  air-air  ejectors,  such  as  considered  in  Section 
XVI,  the  effective  displacement  thickness  is  negative,  i.e., 
the  secondary  flow  is  larger  with  mixing  than  without  mixing 
(Ref.  35).  The  positive  6*  in  the  rocket-air  mixing  system 
can  be  attributed  in  part  to  the  large  heat  release  in  the 
mixing  layer. 


Second-primary  mass  flow  ratio  -  The  experimental  mass 
flow  ratio,  w  /w . ,  is  shown  in  Fig.  26  as  a  function  of 
p^^/p^j .  The  2-D  Core  Theory  predicts  the  correct  form  for 
the  mass  flow  function,  but  the  theoretical  mass  flows  are 
4-5%  higher  than  the  experiment.  The  quantitative  deviation 
between  theory  and  experiment  can  perhaps  be  attributed  to 
deviations  from  conical  flow  in  the  primary  nozzle.  The 
theoretical  computations  were  based  on  a  source-type  flow 
in  the  primary  nozzle ,  but  it  is  well  known  that  the  source 
flow  assumption  can  be  considerably  in  error  for  such  short 
conical  nozzles. 


Mixing  duct  thrust  -  The  mixing  duct  thrust  is  shown  in 

Fig.  27  as  the  ratio  F,/F^  where  F  is  the  vacuum  thrust  of 

the  primary  nozzle.  The  2-D  Core  Theory  predicts  the  thrust 

ratio  fairly  well  at  p„  /p  .  =  0.025,  but  the  theoretical 

*^oa  *^oj 

value  is  about  8%  below  the  experiment  at  ~  0.034. 

Radial  Pitot  pressure  profiles  -  The  radial  distribution 

of  Pitot  pressure  at  the  duct  exit  is  shown  in  Figs.  28a  and 

28b  for  two  values  of  p  „/p  ..  It  should  be  noted  that  the 

*^oa  ^oj 

exit  plane  distributions  are  composite  plots  of  several  short 
duration  tests.  Because  of  the  difficulty  in  exactly  repro¬ 
ducing  all  of  the  experimental  conditions  from  test  to  test, 
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considerable  data  scatter  Is  inevitable  .  The  theory 

adequately  predicts  the  shape  of  the  Pitot  pressure  profile 

for  the  higher  pressure  ratio  case  (Fig.  28a),  but  the 

theoretical  profile  is  displaced  radially  outward  compared 

to  the  experiment.  The  same  is  true  of  the  =  0.025 

case  (Fig.  28b),  and  in  addition,  the  experiment  shows  a 

decrease  in  p  ’  near  the  centerline  of  the  duct.  This  de— 
o 

crease  in  p^'  is  undoubtedly  caused  by  the  presence  of  strong 
shock  waves  in  the  inviscid  core  flow.  The  shock  system 
would  tend  to  be  somewhat  stronger  for  lower  values  of 
p^g^/p^j  compared  to  higher  pressure  ratios  because  the  pri¬ 
mary  flow  is  initially  underexpanded  to  a  greater  extent. 

Radial  composition  profiles  —  Before  considering  the 
experimental  composition  profiles  at  the  duct  exit,  the 
differences  between  the  frozen  and  equilibrium  profiles  will 
be  discussed.  The  experimental  technique  reported  in  Refs.  6  and 
7  was  to  take  a  sample  of  gas,  then  analyze  it  after  removing 
the  water  vapor.  The  composition  profiles  were  then  pre¬ 
sented  as  "dry  mole  fractions"  of  Og,  Hg  and  Ng.  The 
theoretical  dry  mole  fraction  profiles  across  the  mixing  zone 
are  shown  in  Fig.  29.  The  major  characteristic  of  the 
equilibrium  profiles  is  that  the  Ng  tends  to  peak,  and  the 
Hg  goes  to  zero,  near  the  stoichiometric  point  in  the  mixing 
layer.  The  frozen  distributions,  on  the  other  hand,  are 
monotonic  and  continuous  across  the  mixing  zone. 


Measurements  of  Pitot  pressure  and  gas  composition  were 
made  across  the  entire  duct  exit  plane  during  the  AEDC  experi¬ 
ments.  The  flagged  symbols  represent  measurements  made  on 
the  opposite  side  of  the  centerline  from  the  measurements 
represented  by  the  unflu-gged  symbols  (Figs.  28,  30,  34  and  35). 
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Experiinetital  radial  profiles  of  Hg,  O2  Ng  are 

shown  in  Figs.  30a  and  30b  for  p  /p  .  =  0.033  and  0.025, 

oa  oj 

respectively.  Although  there  is  considerable  scatter  in  the 
data,  the  2— D  Core  Theory  with  equilibrium  chemistry  predicts 
the  shape  of  the  experimental  profiles  reasonably  well,  es¬ 
pecially  the  tendency  of  the  Hg  to  go  to  zero  within  the 
mixing  zone.  As  was  the  case  with  the  Pitot  pressure  profiles, 
the  theoretical  profiles  are  displaced  outward  from  the  duct 
centerline  compared  to  the  experimental  profiles.  This  out¬ 
ward  displacement  of  the  mixing  zone  can  be  explained  by  the 
neglect  of  shoch  waves  in  the  inviscid  core  flow.  If  the 
shock  waves  were  included  in  the  theory,  then  the  momentum  of 
the  core  flow  would  he  less  than  the  isentropic  momentum,  and 
the  mixing  zone  would  tend  to  be  oriented  more  toward  the  duct 
centerline.  In  general,  altex-ation  of  the  momentum  ratio  in 
a  two  stream  mixing  process  causes  the  mixing  zone  to  shift 
toward  the  stream  which  is  decreased  in  relative  momentum. 

Comparison  of  1— D  Core  Theory  with  2-D  Core  Theory  -  Al¬ 
though  not  well  adapted  for  analyzing  mixing  systems  operating 
in  the  upstream  choking  mode,  the  1— D  Core  Theory  can  be 
applied  to  such  flows.  The  central  stream  is  assumed  to  be 
Initially  uniform.  Immediately  downstream  from  the  primary 
nozzle  the  primary  and  secondary  streams  are  assumed  to  under¬ 
go  isentropic  area  changes  such  that  the  pressures  are 
equilibrated;  the  combined  primary  and  secondary  stream  areas 
must  equal  the  initial  duct  area.  This  initial  readjustment 
process  between  the  two  streams  is  the  same  as  the  "quasi— 
one— dimensional"  approaches  discussed  in  Section  I,  with  the 
additional  requirement  that  the  readjustment  process  occurs 
instantaneously . 

Applying  the  preceding  procedure  to  the  ducted 
mixing  system,  the  secondary  stream  choking  section  occurs 
at  the  duct  enttance.  For  the  range  of  reported 
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in  Ref.  6,  the  secondary  mass  flow,  w  ,  is  5  to  7%  higher 
than  predicted  by  the  2—D  Core  Theory  (Fig.  26) .  Because  the 
secondary  flow  is  initially  sonic,  =  0.528,  con— 

siderably  lower  than  predicted  by  the  2— D  Core  Theory 
(Fig.  25).  The  wall  pressure  monotonically  decreases  with 
downstream  distance,  and  the  mixing  duct  thrust  is  sub¬ 
stantially  lower  than  predicted  by  the  2— D  Core  Theory. 

These  results  clearly  show  the  need  for  the  2— D  Core  Theory 
lor  mixing  systems  operating  in  the  upstream  choking  mode. 


It  should  be  noted  that  a  less  refined  version  of 
the  1— D  Core  Theory  was  used  in  Ref.  6  to  compute  the  theo¬ 
retical  performance  of  the  rocket— air  mixing  system  operated 
in  the  upstream  choking  mode.  The  flow  upstream  of  the 
secondary  stream  choking  section  was  computed  by  an  approxi¬ 
mate  Inviscid  technique.  The  plume  from  the  primary  nozzle 
was  approximated  by  a  constant  pressure  plume  at  the  secon¬ 
dary  stream  sonic  pressure  (plume  computed  by  a  standard 
method  of  characteristics  program).  Downstream  from  the 
choking  section  the  duct  wall  was  assumed  to  be  shifted 
Inward  to  account  lor  the  difference  between  the  plume  area 
and  the  one-dimensional  primary  stream  area.  The  mixing 
theory  computation  was  then  initiated  at  the  secondary 
stream  choking  section. 


Although  the  preceding  procedure  gives  results 
comparable  to  the  2— D  Core  Theory  for  the  experiments  con¬ 
sidered  here,  it  was  pointed  out  in  Ref.  6  that  the  proce¬ 
dure  is  questionable  for  many  mixing  configurations.  For 
example,  if  the  axial  distance  to  the  secondary  stream 
choking  section  is  large  compared  to  the  duct  length,  then 
most  of  the  influence  of  the  mixing  process  is  neglected. 

2.  Mixing  System  Operated  in  the  Downstream  Choking  Mode 


This  section  presents  the  correlations  of  the  2-D 
Core  Theory  with  experimental  results  for  the  configuration 
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of  Fig.  23b  operated  in  the  downstream  choking  mode  (Ref.  7). 
The  theoretical  computations  for  these  cases  were  made  by 
iterating  on  p^^  until  the  coefficient  determinant,  |  d|  , 
was  equal  to  zero  at  the  duct  exit  plane. 

Wall  pressure  distributions  —  Experimental  wall  pressure 

distributions  are  shown  in  Figs.  31a— c  for  p^„/p^j  of  0.036, 

oa  oj 

0,0304,  and  0.027,  respectively.  The  2— D  Core  Theory  pre¬ 
dicts  the  correct  shape  for  the  pressure  distributions,  but 
the  theoretical  pressure  is  generally  too  high.  This  dis¬ 
agreement  between  theory  and  experiment  is  more  pronounced 

at  p  „/p^  •  =  0.036  than  at  the  lower  pressure  ratios.  As  was 
oa  oj 

the  case  for  the  upstream  choking  mode,  the  deviation  of 
theory  from  experiment  in  the  downstream  part  of  the  flow 
can  be  attributed  to  the  neglect  of  shock  waves  in  the 
theory.  The  inclusion  of  shocks  in  the  theory  would  cause 
the  area  occupied  by  the  inviscid  core  flow  to  increase,  and 
the  area  for  the  subsonic  secondary  flow  would  decrease. 
Consequently,  the  wall  pressures  would  be  decreased. 

Note  that  the  inviscid  core  flow  is  completely 
dissipated  in  this  configuration,  and  the  flow  at  the  duct 
exit  is  in  the  second  mixing  regime.  It  should  also  be  noted 
that  the  flow  is  strongly  influenced  by  viscous  effects. 
Without  mixing,  this  configuration  would  operate  in  the  up¬ 
stream  choking  mode.  Comparing  the  pressure  distributions 
of  Figs.  25  and  31,  one  sees  that  the  addition  of  the  con¬ 
stant  area  extension  to  the  conical  mixing  duct  causes  the 
wall  pressure  to  be  maintained  at  a  high  value  for  the 
entire  duct  length.  This  high  wall  pressure  is  desirable 
for  a  propulsion  system  because  the  wall  pressures  are 
directly  equivalent  to  axial  thrust. 

Secondary— primary  mass  flow  ratio  —  The  experimental  mass 
flow  ratio,  w  /w.,  is  shown  in  Fig.  32  as  a  function  of  the 
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total  pressure  ratio.  The  2— D  Core  Theory  predicts  the 
correct  form  for  the  mass  flow  function,  but  the  theoretical 
mass  flows  are  10—15%  higher  than  experiment.  The  rocket 
engine  used  for  these  experiments  had  a  lower  combustion 
efficiency  than  the  rocket  in  Ref.  6.  Consequently,  the 
theoretical  mass  flow,  w ^ ,  is  about  4%  below  the  experimen¬ 
tal  Wj  for  a  given  •  The  theoretical  and  experimental 

values  of  w  /w .  shown  in  Fig.  32  are  nondimensionalized 
"  J 

with  the  corresponding  value  of  w . .  Because  the  theoretical 

w.  is  lower  than  the  experimental  w.,  the  theoretical  mass 
J  J 

flow  ratio  tends  to  be  too  high.  The  agreement  between 

the  theoretical  and  experimental  values  of  w^  is,  therefore, 

better  by  about  4%  than  indicated  by  Fig.  32. 

Mixing  duct  thrust  —  The  mixing  duet  thrust  ratio, 
is  shown  in  Fig.  33.  Because  the  theoretical  wall  pressures 
are  generally  higher  than  experiment,  the  theoretical  thrust 
ratio  is  also  too  large.  The  theoretical  thrust  ratio  is 
about  6%  larger  than  experiment  at  the  lower  values  of 
p^^/p^ . ,  and  about  16%  too  large  at  the  higher  values  of  the 
total  pressure  ratio. 

Radial  Pitot  pressure  profiles  —  The  exit  plane  distri¬ 
butions  of  the  Pitot  pressure  are  shown  in  Figs.  34a  and  34b 

for  p„-/p_H  of  0,036  and  0.027,  respectively.  Considering 
oa.  oj 

the  data  scatter,  the  agreement  of  theory  and  experiment 
is  fair  for  Pqq^/Pqj  -  0.036.  The  experiment  for  pQ^/Poj  “ 
0.027,  however,  shows  very  low  Pitot  pressures  compared 
with  the  theory  near  the  duct  centerline.  A  strong  shock 
system  has  undoubtedly  occurred  in  the  inviscid  core  flow 
upstream  of  the  duct  exit.  As  was  the  case  for  the  upstream 
choking  mode,  a  stronger  shock  system  would  be  expected  as 
the  total  pressure  ratio  is  decreased  (greater  underexpansion 
of  the  primary  nozzle) . 
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Radial  composition  profiles  —  Exit  plane  profiles  of  H2, 

Oo  and  N„  are  shown  in  Fig.  35  for  p  „/p  .  =  0,036.  The 
z  oa  oj 

agreement  between  theory  and  experiment  is  quite  good  for  this 
case.  The  experimental  Ng  distribution  clearly  shows  the 
peak  predicted  by  the  theory  with  the  assumption  of  equili¬ 
brium  chemistry.  These  results  leave  little  doubt  that  the 
assumption  of  equilibrium  mixing  zone  chemistry  is  good  for 
these  experimental  conditions. 

Exit  plane  choking  parameter  —  These  computations  for  a 
ducted  mixing  system  operating  in  the  downstream  choking 
mode  were  made  by  assuming  that  choking  occurs  at  the  duct 
exit,  and  that  choking  corresponds  to  the  system  of  equations 
becoming  critical  (|d|  =0).  It  is  of  interest  to  compare 
these  results  with  an  exterior  hypothesis  about  the  choking 
of  a  nonuniform  flow.  Herpfer  (Ref.  46)  derived  a  choking 
criterion  for  nonuniform  flow  based  on  the  principle  that 
choking  corresponds  to  a  condition  of  maximum  mass  flux.  He 
assumed  that  there  are  no  radial  pressure  gradients  in  the 
flow  and  that  the  local  spaed  of  propagation  is  the  isen— 
tropic  value.  Choking  then  occurs  when  the  parameter  Ch, 

P  )  v/W^ 

O 

is  equal  to  zero.  This  choking  criterion  states  that  if  a 
portion  of  a  choked  flow  is  supersonic,  then  another  por¬ 
tion  of  the  flow  must  be  subsonic. 

The  parameter  Ch  was  computed  in  the  second 
regime  for  the  ducted  mixing  conditions  presented  in  this 
section.  If  one  defines  an  effective  uniform  stream  Mach 
number,  M^,  such  that 
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C^  =  Tr 


0  -  Muf" ) 

Vp  M  * 


(106) 


then  should  be  unity  when  |d|  =  0  if  Ch  =  0  is  the 
cori’ecl  choking  criterion.  From  the  ducted  mixing  theory, 
was  found  to  be  0.94  ~  0.95  when  |D|  =  0  .  Thus  the 
system  of  equations  lor  the  ducted  mixing  problem  shows  that 
part  of  the  nonuniform  flow  is  indeed  subsonic  when  the  flow 
becomes  critical.  The  value  of  the  secondary  Mach  number, 

M^,  is  even  lower  than  would  be  predicted  by  Ch  =  0.  In¬ 
spection  of  Fig.  31  shows  that  the  theoretical  value  of 
p^/Poa  approximately  0.8  at  the  duct  exit,  therefore, 
is  approximately  0.57. 

3.  Mixing  System  Operated  in  the  Back  Pressure  Dependent  Mode 

This  section  presents  correlation  of  the  2— D  Core 
Theory  with  cxpeirimental  results  for  the  configuration  of 
Fig.  23b  operat'rd  in  the  back  pressure  dependent  mode.  The 
rheoretical  solutions  were  obtained  by  iterating  on  p^^  until 
the  theo’^ct  ical  wall  -piressjre  matched  the  experimental  pres¬ 
sure  ai  the  last  station  In  the  mixing  duct  where  the  static 
pre.ssuri=  was  measured. 


To  svaLaat;;  --q.  (106)  an  effective  value  of  7  must  be 
chosen:  however,  the  r-saits  for  are  insensitive  to  the 
value  chost=n  for  7, 

For  th-,ise  exper Ln’en tal  conditions,  is  higher  than  the 

primary  nozzle  exit  pressure.  The  2— D  Core  Theory  is  not 
applicable  to  an  overexpanded  primary  nozzle,  therefore,  the 
computations  were  riad«  with  the  following  approximation.  The 
conical  primary  nozzle  was  assumed  to  be  shorter  than  the 
experimental  nozzle  (with  a  smaller  area  ratio)  such  that  the 
nozzle  exit  pressure  was  slightly  higher  than  ' 
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Wall  pressure  distributions  -  Axial  distributions  of 

wall  static  pressure  are  shown  in  Figs.  36a  and  36b  for  two 

values  of  p  /p  . .  In  these  cases  the  theoretical  pressure 
^oa  *^oj 

naturally  agrees  well  with  experiment  in  the  downstream  part 
of  the  duct,  but  the  theory  is  too  low  near  the  duct  entrance. 
As  a  result  of  the  low  theoretical  values  for  the 

theoretical  values  of  the  secondary— primary  mass  flow  ratio 
are  too  high  (by  32%  and  28%,  respectively).  The  theoreti¬ 
cal  thrust  ratios  are  much  less  in  error,  being  low  by  6—7%. 
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SECTION  XVIII 
CONCLUSIONS 

An  analytical  model  has  been  presented  for  the 
turbulent  mixing  of  coaxial  streams  inside  a  duct  of  arbi¬ 
trary  shape.  The  influences  of  nonuniform  inviscid  core  flow 
and  of  equilibrium  chemical  reactions  have  been  included  in 
the  analysis.  Even  though  the  current  detailed  knowledge 
about  turbulent  flows  with  chemical  reactions  is  meager,  the 
use  of  integral  methods  permits  reasonably  accurate  compu¬ 
tations  of  the  flow  in  complex  mixing  systems.  An  important 
point  in  the  analytical  model  is  the  technique  used  to  couple 
the  integral  solution  for  the  mixing  layer  and  secondary  flow 
with  the  solution  for  the  nonuniform  inviscid  core  flow 
(treated  with  the  method  of  characteristics) .  This  technique 
should  be  applicable  to  other  supersonic  viscous  interaction 
problems  in  which  integral  methods  can  be  applied  to  the 
viscous  layer. 

The  analytical  and  experimental  results  which  were 
discussed  in  Section  XVII  clearly  show  that  viscous  effects 
are  predominant  in  certain  types  of  ducted  mixing  systems. 

In  some  cases,  the  viscous  effects  completely  change  the 
character  of  the  flow  from  that  predicted  by  inviscid  theories. 
Consequently,  it  is  essential  that  the  viscous  and  inviscid 
portions  of  the  flow  be  computed  simultaneously  for  such 
mixing  systems. 

The  main  emphasis  of  this  investigation  has  been 
placed  on  predicting  the  performance  of  ducted  mixing  systems 
which  are  strongly  influenced  by  viscous  effects.  Of  course, 
the  theoretical  model  is  also  applicable  to  flows  which  are 
weakly  influenced  by  mixing,  such  as  supersonic  air— air 
ejectors  operating  in  the  upstream  choking  mode.  The  present 
theory  accurately  predicts  the  experimental  performance  of 
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such  air— air  ejectors,  but  the  theory  offers  no  quantitative 
improvement  over  the  superposition  technique  of  Chow  and 
Addy  (Ref.  5). 

Many  assumptions  have  been  made  in  the  development 
of  the  present  analytical  model.  Extensive  correlations  of 
the  theory  with  experiments,  however,  have  shown  that  the 
various  assumptions  are  acceptable,  or  at  least  not  pro¬ 
hibitively  inaccurate.  The  method  of  treating  the  flow  in 
the  mixing  layer  appears  to  be  satisfactory  considering  the 
current  state  of  knowledge  about  free  turbulent  mixing  in 
the  regime  of  interest.  There  is  no  point  in  discarding  the 
assumptions  of  mixing  zone  profile  similarity  and  of  unity 
Prandtl  and  Lewis  numbers  until  detailed  and  precise  experi¬ 
mental  results  are  available  for  mixing  with  large  density 
gradients  and  chemical  reactions.  The  assumption  of  equili¬ 
brium  mixing  zone  chemistry  yields  good  results  for  the 
rocket— air  mixing  systems  considered  in  this  investigation. 
However,  the  theory  should  be  applied  with  caution  to  mixing 
systems  for  which  there  is  no  experimental  verification  for 
the  assumption  of  equilibrium  chemistry. 

The  weakest  aspects  of  the  present  theory  are  (1) 
the  model  for  the  eddy  viscosity,  and  (2)  the  assumption  that 
the  inviscid  core  flow  is  irrotational  (without  shocks) .  The 
eddy  viscosity  model  is  based  on  the  Prandtl  incompressible 
model,  which  is  known  to  be  inadequate  when  the  secondary 
velocity  is  larger  than  about  0.3  times  the  primary  velocity. 
Although  the  correction  used  to  account  for  the  influence 
of  density  gradients  gives  satisfactory  results  for  the  cases 
considered  in  this  investigation,  the  correction  needs  to  be 
checked  against  a  broad  range  of  experiments.  It  is  hoped 
that  the  current  research  on  free  turbulent  mixing  will  re¬ 
sult  in  an  improved  semi— empirical  model  for  the  eddy 
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viscosity;  such  an  improved  eddy  viscosity  model  could 
easily  be  incorporated  into  the  present  ducted  mixing  theory. 

It  is  thought  that  the  major  source  of  quanti¬ 
tative  disagreement  between  the  present  theory  and  ducted 
rocket— air  mixing  experiments  is  the  neglect  of  shock  waves 
in  the  inviscid  core  flow.  There  is  no  conceptual  problem 
in  incorporating  the  rotational  method  of  characteristics 
(with  shock  waves)  into  the  present  theory.  The  method  of 
computing  the  inviscid  core  flow  (Section  XIII)  would  be 
the  only  change  required  in  the  existing  computer  program; 
the  system  of  differential— integral  equations  for  the  mixing 
layer  and  the  inviscid  secondary  flow  would  be  unchanged. 

In  spite  of  the  deficiencies  in  the  theoretical 
model  developed  during  this  investigation,  the  results  are 
certainly  qualitatively  correct.  In  addition,  the  theoreti¬ 
cal  results  are  sufficiently  accurate  so  that  the  theory 
may  be  considered  a  useful  tool  for  the  engineering  analysis 
of  ducted  mixing  systems. 
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APPENDIX  I 

AUXILIARY  FUNCTIONS  FOR  1-D  CORE  THEORY 

In  Section  VII,  the  basic  differential-integral 
equations  for  the  ducted  mixing  problem  were  transformed 
into  a  system  of  linear  first  order  differential  equations 
which  is  amenable  to  numerical  solution.  The  numerous 
auxiliary  equations  needed  to  evaluate  the  coefficients 
F,  G  and  H  are  presented  in  this  appendix. 

The  following  equations  are  developed  in  logical 
order  in  that  all  the  parameters  in  a  particular  equation 
are  developed  prior  to  use  in  that  equation.  All  of  the 
equations  for  the  parameters  which  appear  directly  in  the 
coefficients  F,  G  and  H  are  denoted  by  an  asterisk  at  the 
equation  number,  for  example  eq.  (A— 11)  . 

Note  that  in  the  evaluation  of  the  terms  for  the 
coefficients  the  independent  variables  are  considered  to  be 
(x,  p^,  r^,  b)  in  the  first  regime,  (x,  p^,  u^,,  b)  in  the 
second  regime  and  (x,  p  ,  u  ,  u  )  in  the  third  regime. 

Vf  C  W 

1.  Isentropic  Stream  Equations 

A  number  of  relations  which  are  used  in  all  three 
regimes  are  obtained  from  the  assumption  that  the  Inviscid 
streams  are  isentropic.  The  following  parameters  are 
assumed  to  be  prescribed:  p^^ ,  T^^ ,  R^,  Rj ,  Cp^, 

c  . . 

PJ 

From  the  definition  of  the  ratio  of  specific  heats: 


V 


CpA.. 

Cpa,“ 


(A-1) 


(A-2) 
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From  the  isentroplc  pressure  ratio  equation: 


Differentiating  eqs.  (A-3)  and  (A-4)  with  respect  to  p^: 


1 

Mol* 

CL 

1 

ri-i 

d  Pw 

Pw  Mj 

From  the  energy  equation  : 


(A-3) 


(A— 4) 


(A-5) 


(A-.6) 


(A-7) 


Use  of  the  energy  equation  (and  eqs.  A— 3  and  A-4)  in  this 
form  implies  that  the  specific  heats  are  invariant  over  the 
range  of  temperature  from  0  to  Tq,  whereas  in  Section  V  the 
total  temperature  was  defined  by  assuming  that  Cp  is  constant 
from  T  to  Tp.  For  rocket— air  mixing,  no  error  is  involved 
in  the  secondary  stream  which  is  relatively  low  terapera+ure. 
For  the  rocket— air  mixing  experiments  of >Refs.  6  and  7, 

Tj  1100“K,  and  the  error  in  Tj  caused  by  using  eq.  (A— 8) 
is  less  than  1%. 
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(A-8) 


Differentiating  eqs.  (A-7)  and  (A-8)  with  respect  to  p^: 


(j—  Ma/TX* 

d  M<y 

d  Pvo  ” 

apvi 

JKAj 

dp,^ 

d  pw 

From  the  definition  of  the  Mach  number; 

1/ 

Uj  =  (rj  Ti)  ^ 

Differentiating  eqs.  (A-11)  and  (A-12)  with  respect  to  p„: 

w 


(A-9) 


(A-10) 


(A-ID* 


(A-12)* 


d  Moa  dpvg  2-4^^  ^  ^dpw 


d 

d  p*j 


dPw 


+  4-  M  j  rj  ^ j 


«.i 


dTj 

dpw 


CA-14)* 
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From  the  perfect  gas  law: 


Differentiating  eqs.  (A— 15)  and  (A— 16)  with  respect  to  p  : 


(A-15)* 

(A-16)* 


df«. 

- 

ciTa, 

1 

1 

Pw  ~ 

To. 

dPw  ^ 

_ 

1 

= 

Tj 

-  4“ 

d 

2.  Duet  Wall  Equations 


(A-17)* 


<A-18)* 


The  duct  wall  Is  assumed  to  be  represented  by  the 
following  pair  of  equations,  where  the  constantsCj^  -  and 
are  specified: 

'Tn  =  C|  -I-  X-  +  C^.X‘  +  C5X, 

pc  ^  TCq^  (A-19a) 


(A-19b) * 


Differentiating  eqs.  (A-19a)  and  (A-19b)  with  respect  to  x: 


d  TVi 

_ V 


=  “h  Z  ^3  ZC  4-  4-  ^  ^  pC 

(A-20a) * 


<L7C 


7 


(A-20b)* 
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3.  Mixing  Zone  Property  Equations 

The  mixing  zone  total  temperatures,  specific  heats 
and  gas  constants  are  prescribed  (Section  V): 


To  To  CC) 

(A-21) 

Cp  =  CpCc) 

(A-22) 

R  =  RCS) 

(A-23) 

These  equations  are  also  differentiated  with  respect  to  C 
to  provide  dT  /dC,  dc  /dC  and  dR/dC. 

4.  First  Regime  Equations 

From  Section  V,  the  equations  for  the  mixing  zone 
velocity  and  composition  are: 


^  =.  COS  IT  — ^  ^  ^  (A_2S)* 


Differentiating  eqs.  (A-24)  and  (A~25): 

^  =o 


SPw  b  Adpw  dpgjJ 


<^a 


(A-26) 


(A-27) 


(A-28) 
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(A-29)* 


(A-30) 


(A-31) 


(A-32) 


(A-33) 


From  eqs.  (A-21)  -  (A-23)’one  obtains  dT^/dC,  dCp/dC,  and 
dR/dC ,  and 


^“To 

hCf  _  - 

(A-34) 

£[o  ^ 

Cp  ^  ^ 

CA-35) 

Pw  ^Pvg 

^To  _ 

<i  C 

(A-36) 

^C.p 

SC.  Cp 

dC. 

(A-37) 

= 

^C.  d  ^ 

dC 

(A-38) 
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C^To 

dTo 

as 

^  Cp  ^ 

1 

a  c.f» 

^  b 

dCT 

_ 

. 

Al. 

^  b 

dC 

The  mixing  zone  temperature,  T,  is  given  by  the  energy 
equation: 


Differentiating  eq.  (A-42) ,  one  obtains 


O 


2.0^  \ 

_  UL 

^  P\M 


'ZLLC^ 


(A-39) 


(A-40) 


(A-41) 


(A-42) 


(A-43) 


(A-44) 


(A-45) 
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‘SEy 


(A-<46) 


The  mixing  zone  density  Is  given  by 

p  =  Pw/rT 

Differentiating  eq.  (A-47) ,  one  obtains 


CPw 


(A-47)* 

(A-48)* 

(A-49)* 
(A-50) * 
(A-51)* 


Equations  (A— 24)  —  (A— 51)  completely  specify  the  terms 
which  appear  in  the  mixing  zone  integrals  of  the  coeffi¬ 
cients  F,  G  and  H. 

The  mixing  zone  shear  stress,  is  given  by 

eq.  (44a)  of  Section  VI; 


(44a)* 
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where  the  half-radius  density,  Is  evaluated  from  eq.  (A-47) 
at  C  =  J  and  i  (uj  +  ‘  first  regime  mixing  zone 

constant,  k^,  is  given  by  the  modified  Donaldson  and  Gray 
equation  (eq.  42,  Section  VI), 


s  a. 007(0.  £4  -h  0. 34-  eif  4g.  M^^)j 


(A-52) 


where  M  is  the  Mach  number  corresponding  to  u  and  C  = 

m  "  m  ■* 

In  calculating  M^,  the  variables  T  ,  and  c^«  are  deter- 

m  m  m  pm 

mined  by  eqs.  (A— 21)  —  (A-23)  and  (A-42) ,  The  definition 

for  7^  (7_  =  c _ /(c _  —  R_))then  permits  determination  of 

m  .m  pm  pm  m  ' 

from: 

m 


^  w  ~  C^fn  ^iri 


5,  Second  Regime  Equations 

The  mixing  zone  velocity  is  given  by  (SectionyV) : 

CC  —  +  C  0  S  “TT 


(A-53) 


(A-54) 


The  mixing  zone  composition,  C,  is  given  by 


C  ^Cc.  Cl  4  0,0 


(A-55) 


where  C  is  determined  by  eq.  (25a)  (Section  V) : 
c 


LLj^Ucu 


(A-56) 


The  velocities  u  and  u.  are  given  by  eqs.  (A— 11)  and 

“  J 

(A-12) .  Differentiating  eq.  (A— 54): 


^IL 


=  o 


(A-57) 
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1 


-h  C.OS  77~  ^  J 


£\ 

cos  TT  / 


dU-a^ 

d  P« 


Differentiating  eq.  (A— 56) ,  one  obtains 

^  =o 

hX- 

dCc  _  _J _ 

<^Pw  «-j) 

dCc  _  ^ 

Differentiating  eq.  (A-55),  one  obtains 


o 


iPw 


CoS  TT  b) 


^Cc 


(A-58)* 


(A-59) * 


(A-60) * 


(A-61) 


(A~62) 


(A-63) 


(A-64) 


(A-65) 
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Equations  (A-21)  -  (A-23),  (A-34)  and  CA-39) 
(A*^l)  are  applicable  in  this  regime.  In  addition, 

dTo 

^  Pw  ^Pvo  ^  5 

gl 

^  Pnn  ci  C 


i  psV  S  Pvj  ci  2 

d 

^ do 

^  U.C.  ^  d  C 


(A-66) 


(A-67) 


(A-68) 


(A-69) 


(A-70) 


(A-71) 


(A-72) 


(A-73) 
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The  mixing  zone  temperature  is  given  by  eq.  (A-42) . 
Differentiating  eq.  (A— 42) ,  one  obtains  eq.  (A— 46)  and 


^  =o 

nr  J—  /^,je  ^  i£PL\ 


(A-.74) 


(A-75) 


^  ^C. 


LL^ 

Uc/ 


(A-76) 


The  mixing  zone  density  is  given  by  eq.  (A-47) . 
Differentiating  eq.  (A— 47) ,  one  obtains  eqs.  (A— 48)  and 
(A-51)  and 


—  =  +  T  1^') 


(A-77)* 


iL  4- 

^ac. 

The  mixing  zone  shear  stress  is  given  by  eq.  (44b) 
of  Section  VI: 


(44b)* 


The  density,  is  calculated  at  u^  ■=  i(u^  +  u^>  and 

-  1/2  . 
m  c 

The  second  regime  mixing  constant,  ^hj  given  by 
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O,  0//[o.£6-hAS4-sxf>^3.4£ 


(A-79) 


Where  H  Is  calculated  with  eq.  (A— 53) . 
in 

6.  Third  Regime  Equations 

The  mixing  zone  velocity  is  given  by  (Section  V) : 


=  'h  ^  COS  TT 


CA-80) 


The  mixing  zone  composition  is 


C  =-Cvg  !^C^c-C^)(l+  cesTT^') 


(A^Sl) 


where 


— — 

C.  ^ 


and 


U-J  -  Ua^ 


C-  ^ 


CL  *  ft-r 

a/ 


(A-56) 


(A*>82) 


The  velocities  u  and  u.  are  given  by  eqs,  (A— 11)  and  (A— 12) . 

a  j 

Differentiating  eq.  (A— 80) : 


cA<V/ 


A7£. 

(A-83)* 


^  iXi 


(A-84)* 
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Differentiating  eq.  (A-56) ,  one  obtains 

^C.e.  _  _  Q 

^  U.M 

^Cc  *  [V  ^\dLLtu  —  di4j' 

_  _| _ 

^  U-  c.  C^j  "*■ 

Differentiating  eq.  (A--82) ,  one  obtains 

^  ^C.vj  __  Q 

X  d  U<C 


__j _ 

Cua.-  O 


^Cw  ^  -J _ 

^  U.W 


(A-85)* 


(A-86)* 


(A-87) 


(A-88) 


(A-89) 


(A-.90) 


(A-91) 


(A-92) 
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Differentiating  eq.  (A-81) ,  one  obtains 


dc 


±(j+cesTr^(^ 


^gc 


<iPW/ 


iSUt  ^  ' 


^5e 

^  U»c 


<^c> 

^U.w 


^(l-cesTT^')!^ 

'<^Uw 


Equations  (A-21)  -  (A— 23)  are  applicable  to  the 
third  regime.  The  derivatives  of  the  mixing  zone  gas 
properties  are: 


^  aTo 

d  C. 

Pw  d  C 


im  =.  ^  dT^o 
^x  “  cix  at 

^  dCp 

^P>M  ^Pw  Jc* 


(A-93) 

(A-94) 

(A-95) 

(A-96) 

(A-97) 

(A-98) 

(A-99) 

(A-lOO) 

(A-101) 
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^Cp  _  ^  C.  d  cp 

U-c-  d  c 

(A-102) 

^  d  Cp 

(A-103) 

^  I?  —  df? 

dC. 

(A-104) 

—  bC.  d 
^  Pw  d  C 

(A-105) 

dR 

^U-c  ^^CdCl 

(A-106) 

42.  ,  ^c.  dR 
^  U.WJ  ^  U.VW  d  C 

(A-107) 

The  mixing  zone  temperature  is  given  by  the 


energy  equation  (eq.  A— 42) ♦  Differentiating  eq.  (A— 42)  one 
obtains  eqs,  (A-75)  and  (A-76)  and 


I 

zc 


74 


^  16 
^7C 


(A-108) 


^Te  ^  ^  fj  r  ^  ^  IJ  ^  (A-109) 

2.Cp^\  P  dU-vJ 

The  mixing  zone  density  is  given  by  eq.  (A— 47) . 
Differentiating  eq.  (A-47) ,  one  obtains  eqs.  (A-77)  and 
(A— 78)  and 


_ /i?  il  4- T 


(A-110)* 
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_ 4,-r 

diAwJ 


(A-IID* 


The  mixing  zone  shear  stress  Is  given  by  eq.  (44c) 
of  Section  VI, 

^  (44c)* 

The  density,  p  ,  is  calculated  at  u  =  J(u  +  u  )  and  C  = 

_  _  "^m  Di  c  w  m 

i(C^  +  C^) .  The  third  regime  mixing  constant,  kjjj,  is  the 
same  as  and  is  calculated  with  eq.  (A— 79) . 
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APPENDIX  II 

AUXILIARY  FUNCTIONS  FOR  2-D  CORE  THEORY 


'  In  this  appendix  the  auxiliary  functions  required 

for  evaluation  of  the  coefficients  F,  G  and  H  (Section  XII) 
are  presented.  Host  of  the  auxiliary  functions  are  the  same 
as  those  for  the  1— D  Core  Theory  (Appendix  I).  Consequently, 
only  those  functions  which  differ  from  the  1-D  Core  Theory 
are  presented  here.  As  in  Appendix  I,  the  equations  for  the 
parameters  which  appear  directly  in  the  coefficients  F,  G 
and  H  are  denoted  by  an  asterisk  at  the  equation  number . 

As  was  the  case  in  evaluating  the  auxiliary 
functions  for  the  first  regime  of  the  1— D  Core  Theory,  the 
variables  (x,  p^,  r^,  b)  are  considered  to  be  independent. 

1.  U«ntrepic  Stream  Equations 

The  secondary  stream  velocity,  u. ,  is  calculated 

A 

from  eq.  (A-11) .  The  invlscid  core  velocity  at  r^,  V^, 

is  computed  from  eq.  (A— 12)  by  replacing  u.  by  V.. 

*J  ^ 

Similarly,  dV^/dp  is  calculated  with  eq.  {A-14) .  The 
density  at  r^,  is  obtained  from  eq.  (A— 16)  by  replacing 

by  p . .  Similarly,  dp. /dp  is  calculated  from  eq.  (A— 18) . 

J  1  1  W 

2.  Velocity  Profile  Equation 

The  mixing  zone  velocity  is  calculated  from 
((is  (A-24a)* 

where 

—  V^'  Cos*  (A-112)* 
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The  normal  velocity  component  at  is  given  by 

=■  V:  s  ‘  «  ®t 

Differentiating  eq-  (A— 24a) ,  one  obtains 


CCS  TT 


b  /  s  56 


(A-113)* 

<A-114)* 


^  -k0+  oes  TT 


A&p« 


—  I 

dpvj/  c^FW 


(A-115)* 


Equations  (A— 114)  —  (A— 116)  contain  partial  derivatives 
Su./Sx,  Su./Sp  ,  and  Su./Sr.  which  must  be  obtained  from 

1  X  iT  1  X 

the  method  of  characteristics  solution  of  the  inviscid 
core  flow  (Section  XIII).  These  derivatives  are  presented 
in  the  next  section. 

3.  Partial  Darivatives  at  Inner  Mixing  Zone  Boundary 

One  obtains  from  the  method  of  characteristics 
solution  (Section  XIII)  the  partial  derivatives  3u/3x, 
3u/dr,  5p/dxj  and  dp/Sr.  It  is  important  to  recognize  that 
these  derivatives  are  obtained  from  a  system  in  which  the 
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independent  variables  are  x  and  r.  Thus,  these  derivatives 

will  be  denoted  by  1  1  with 

the  subscript  denoting  the  variable  held  constant  during 
the  differentiation. 

The  derivatives  appearing  in  eqs.  (A— 114)  —  (A— 116) 
are  used  in  a  system  where  the  independent  variables  are 
Cx,  p^,  r^),  thus 


^ILt 


One  can  write,  in  general 


^  1>L: 


and 


1  —  ^ 

.  ^ 

lx*' 

k.r 

ar 

-  ^ 

+•  ^ 

iE. 

r  6X 

'fjp 

ax 

Thus  one  can  compute  the  necessary  derivatives  from 


^u. 

r  apw 

x,ri  dx 

and 


au.; 

^1 

iE.) 

ar- 

Pw 

ar  lx” 

x.,r-  arl 

The  derivative 


.must  now  be  evaluated  to 


(A- 


(A- 


solve  eqs.  (A— 118)  and  (A-119) .  Differentiating  eq.  (A-112) : 


118) 


119) 
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=,  -V/  stn  ©/  ^ — 

;£,rc  ‘  ^Pw 


TC^^i  oPv/ 


(A-120) 

m 


From  the  Bernoulli  equation, 


^Pwg 


I 


_  Ay  I  _ 


(A-120a) 


The  derivative  is  given  by 

^  _  >1^1*’-  I 

(iPw  a:,ri  P'V'-'®' 

The  results  of  eqs.  (A-118)  —  (A— 120)  allow  the  solution  of 
eqs.  (A-114)  -  (A-116)  for  the  necessary  derivatives. 


(A-120b) 


4.  Mixing  Zon«  T«mperature 

The  mixing  zone  temperature  is  assumed  to  be  given 

by  eq.  (A— 42) .  To  be  correct,  the  kinetic  energy  should  be 

2  2 

represented  by  V  rather  than  u  ;  therefore,  use  of  eq. 
(A-42)  implies  that  the  flow  angles  in  the  mixing  zone  are 
not  large,  i.e.,  cosS^i;  1.  As  discussed  in  Section  X,  the 
system  of  differential— integral  equations  was  derived  on 
the  basis  that  cos  6  ii:  1;  therefore,  use  of  eq.  (A— 42)  does 
not  further  limit  the  validity  of  the  2— D  Core  Theory. 

All  of  the  derivatives  of  the  temperature  are 
the  same  as  in  Appendix  I,  with  the  following  exception; 


^u. 


CA-121) 
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5.  Mixing  Zone  Shear  Stress 

The  shear  stress,  Is  given  by 

(A-122) 

where  Is  evaluated  at  u  ,  and 

+  (A-123) 
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APPENDIX  III 

INPUT  AND  PRINTOUT  PARAMETERS  FOR  COMPUTER  PROGRAMS 


This  appendix  lists  the  parameters  which  must  be 
prescribed  for  operation  of  both  the  1— D  Core  Theory  and 
2-D  Core  Theory  computer  programs.  The  flow  field  variables 
which  are  printed  out  are  also  listed. 

1.  1-D  Core  Theory  Program 

Input  parameters  —  The  following  parameters  must 
be  prescribed: 

1.  Incompressible  mixing  zone  constants, 

■'oil’  ''oiir 

2.  Integration  step  size,  Ax,  and  maximum  x  (L) . 

3.  Primary  nozzle  radius,  r^. 

4.  Initial  static  pressure, 

5.  Constants  in  duct  wall  equation,  C^^  -  C^,  x^. 

6.  Inviscid  stream  stagnation  conditions, 

^’oj  ’  '^oa  ’  "^o  j  ■ 

7.  Inviscid  stream  gas  properties,  c  c  R  , 

pa  pj  a 


8.  Mixing  zone  property  functions;  T^^CC) ,  Cp(C), 

R(C)  (an  array  of  20  values  from  C  =  0  to  C  =  1) . 

Printout  parameters  —  The  following  parameters  are 
printed  out  at  each  axial  station  in  the  first  regime 
computation: 

1.  Wall  pressure,  p^. 

2.  Secondary  stream  Mach  number, 
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3.  Primary  stream  Mach  number,  Mj . 

4.  Width  of  mixing  zone,  b. 

5.  Radius  of  inner  mixing  zone  boundary,  r^. 

6.  Mixing  zone  growth  rate,  db/dx. 

7.  Total  stream  thrust. 

8.  (r^  +  b)/r^. 

9.  Species  conservation  parameter,  Q. 

10.  Coefficient  determinant,  |  d]  . 

The  following  parameters  are  printed  out  at  each 
axial  station  In  the  second  regime  computation: 

1.  Wall  pressure,  p  . 

W 

2.  Secondary  stream  Mach  number,  M^. 

3.  Width  of  mixing  zone,  b. 

4.  Growth  rate  of  mixing  zone,  db/dx. 

5.  Centerline  velocity  ratio,  u  /u  . . 

C  J 

6.  Centerline  composition,  C  . 

7.  b/r^. 

8.  Total  stream  thrust. 

9.  Species  conservation  parameter,  Q. 

10.  Coefficient  determinant,  |d| . 

11.  Choking  parameter,  Ch. 

The  following  parameters  are  printed  out  at  each 
axial  station  in  the  third  regime  computation: 

1.  Wall  pressure,  p^. 

2.  Centerline  velocity  ratio,  u  /u . . 

C  J 


185 


AEDC-TR-6B-270 


3.  Centerline  composition,  C^. 

4.  Wall  velocity  ratio, 

5.  Wall  composition,  C^- 

6.  Centerline  Uach  number, 

c 

7.  Wall  Mach  number,  M^. 

8.  Total  stream  thrust. 

9.  Species  conservation  parameter,  Q. 

10.  Coefficient  determinant,  |d|  . 

11.  Choking  parameter,  Ch. 

2.  2-D  Core  Theory  Program 

Input  parameters  —  The  input  parameters  are  the 
same  as  for  the  1— D  Core  Theory  program,  with  the  addition 
that  the  flow  conditions  across  the  primary  nozzle  exit 
plane  must  be  prescribed.  An  array  (typically  11  points)  of 
(p/P^jj  0)  vs.  r  is  used  (for  0  <  r  <  r^^) . 

Printout  parameters  —  The  following  parameters  are 
printed  out  at  each  axial  station  in  the  first  regime  compu¬ 
tation: 

1.  Wall  pressure,  p  , 

w 

2.  Secondary  stream  Mach  number,  M^. 

3.  Mach  number  at  inner  boundary  of  mixing 
zone, 

4-  Width  of  mixing  zone,  b. 

5.  Radius  of  inner  mixing  zone  boundary,  r^. 

6.  Mixing  zone  growth  rate,  db/dx. 

7.  Total  stream  thrust, 

8.  (r^  +  b)/r^. 
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9.  Species  conservation  parameter,  Q. 

10.  Coefficient  determinant,  |d|  . 

11.  Flow  angle  at  inner  mixing  zone  boundary,  0^. 

One  also  has  the  option  of  printing  out  a  radial  array  of 
(p/p^j ,  9)  across  the  inviscid  core  flow. 

The  printout  parameters  in  the  second  and  third 
mixing  regimes  are  the  same  as  for  the  1—D  Core  Theory 
program. 
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13  *BSTR*cT^  theoretical  investigation  of  the  ducted  turbulent  mixing  of  co 
axial  streams  is  presented.  The  effects  of  equilibrium  chemical  reac¬ 
tions  in  the  mixing  layer  are  included  in  the  analysis.  In  Part  I,  the 
theory  is  developed  for  ducted  mixing  systems  in  which  the  inviscid 
portions  of  the  flow  are  assumed  to  be  one-dimensional .  The  integral 
method  is  used,  and  the  mixing  zone  velocity  profiles  are  assumed  to  be 
similar.  Turbulent  Prandtl  and  Lewis  numbers  are  assumed  to  be  unity, 
and  the  turbulent  shear  stress  along  one  control  surface  in  the  mixing 
layer  is  computed  by  use  of  a  semi-empirical  model  for  the  eddy  viscos¬ 
ity.  The  theory  is  formulated  to  include  the  downstream  regime  where 
the  mixing  layer  extends  across  the  entire  duct.  The  system  of  inte¬ 
gral  equations  is  transformed  into  a  form  suitable  for  numerical  solu¬ 
tion  by  standard  methods.  In  Part  II,  the  theory  is  extended  to  in¬ 
clude  the  influence  of  a  nonuniform  supersonic  Inviscid  core  flow.  A 
technique  is  developed  for  coupling  the  numerical  solution  of  the  in¬ 
viscid  core  flow  (from  the  method  of  characteristics)  with  the  solutio 
of  the  integral  equations  for  the  mixing  layer.  The  theory  is  correla 
ted  with  low  speed  experiments,  with  supersonic  air-air  ejector  experl 
ments  and  with  experiments  on  ducted  rocket-air  mixing.  The  correla¬ 
tions  show  that  the  theory  gives  satisfactory  results  for  the  perform¬ 
ance  of  ducted  rocket-air  mixing  systems  which  are  representative  of 
air  augmented  rockets. 
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